load packages:
library(googlesheets4) #scrapes google sheets
library(dplyr) # dataframe manipulation
library(ggplot2) # plotting package
library(lubridate)
library(tidyr)
library(readr)
library(scales)
library(growthrates)
library(patchwork)
library(grid)
library(qiime2R)
library(phyloseq)
library(viridis)
library(magrittr)
library(vegan)
library(stringr)
library(ape)Create a color pallete
pp_palette <- c(
# Original + Extended Palette (50 colors)
"#7D8F69", "cornflowerblue", "#F5C45E", "#C1D3D5", "plum",
"hotpink4", "lightsalmon", "steelblue", "indianred", "springgreen4",
"#102E50", "maroon", "lightgreen", "lightblue", "gold",
"#BFA5A0", "#A26769", "blue", "red", "lightpink",
"#71816D", "#FFE6E6", "#E5C1CD", "#A3B18A", "salmon", "magenta", "orchid","deepskyblue", "tomato",
"cyan", "hotpink", "yellow", "orange",
"#C0A98E", "#7A9E9F", "#A1869E", "#BFD8B8", "#F7DD72","green", "darkgreen","purple",
"#FFE5AD",
"#CDB4DB"
)
scales::show_col(pp_palette)Import sequence data from QIIME 2: Amplicons were processed on a QIIME 2 pipeline, which includes DADA denoising. Seqeunces were classified with Native Bayers Classifier
Import and prepare metadata:
meta <- read.csv("SampleMetaData.csv")
meta$SampleID <- paste0(meta$SampleID, "_6283")
row.names(meta)<- meta$SampleID
meta$TreatRep <- paste(meta$PCRDate_Treat, meta$SampleInfo_Rep, sep = "_")
#keep only the P. bahamense project
pyro_only <- meta[grepl("^Pyro", meta$Project, ignore.case = TRUE), ]
META <- sample_data(pyro_only)Import and prepare taxonomy data:
taxonomy <- read.csv("taxonomy.csv", stringsAsFactors = FALSE)
names(taxonomy) <- c("row", "tax", "Confidence")
row.names(taxonomy) <-taxonomy[[1]] #move the feature ID column to become the row names
taxonomy <- taxonomy[,(-1)] #delete the feature ID column
taxonomy <- separate(taxonomy, tax, c("Domain","Phylum", "Class", "Order", "Family", "Genus", "Species", "D7", "D8", "D9", "D10", "D11", "D12", "D13", "D14"), sep = ";", fill = "right")
taxonomy <- taxonomy[,c(1:7)] #keep first seven levels
taxonomy$D0 <- with(taxonomy, ifelse(Order == " o__Chloroplast", "Chloroplast", "Bacteria")) #create column for ordering
# Reorder columns and turn into phyloseq tax_table
col_order<- c("D0", "Domain","Phylum", "Class", "Order", "Family", "Genus", "Species")
taxonomy<- taxonomy[, col_order]
taxmat <- as.matrix(taxonomy)
TAX = tax_table(taxmat)Subsetting the phyloseq object:
#Keep only Pyro samples
ps <- prune_samples(rownames(pyro_only), ps)
#Remove OTUs that have zero counts in all Pyro samples
ps <- prune_taxa(taxa_sums(ps) > 0, ps)
ps <- merge_phyloseq(ps, TAX, META)Bar plot at domain level:
#Keep only Bacteria + Chloroplast
ps<- subset_taxa(ps, D0 %in% c("Bacteria", "Chloroplast") )
#Convert counts to relative abundance (%)
psra<- transform_sample_counts(ps, function(OTU) 100* OTU/sum(OTU))
#Group OTU by DO level
glomD1<- tax_glom(psra, 'D0')
#plot RA
taxabarplot <- plot_bar(glomD1, x = "SampleID", fill = "D0") +
scale_y_continuous(expand = c(0, 0)) +
theme_classic() +
theme(
legend.title = element_blank(),
text = element_text(size = 14),
axis.text.x = element_text(angle = 90)
) +
ylab("Relative Abundance (%)") +
xlab("") +
scale_fill_manual(values = c("lightgrey", "#7BB03B")) +
facet_grid(~Project, scales = "free_x")
taxabarplotRemove chloroplast sequences from data:
Calculate observed richness:
plugin <- ps_nochloro %>%
estimate_richness(measures = "Observed") %$% Observed
#Pull project metadata field from phyloseq object
Project <- ps_nochloro %>% sample_data %$% Project
#Combine richness + project into one dataframe
richness<- data.frame(plugin, Project)
names(richness) <- c("richness", "Project")
richness %>%group_by(Project) %>% summarize(mean = mean(richness), min = min(richness), max = max(richness))## # A tibble: 1 × 4
## Project mean min max
## <chr> <dbl> <dbl> <dbl>
## 1 Pyro 97.5 57 295
Relative abundance and sums the counts of all OTUs that belong to the same Order:
ps_nochloro_RA<- transform_sample_counts(ps_nochloro, function(OTU) 100* OTU/sum(OTU))
ps_nochloro_RA_glomO<- tax_glom(ps_nochloro_RA, 'Order')Bar plot at the order level:
# subset by pyro samples
pyro <- subset_samples(ps_nochloro_RA_glomO, Project == "Pyro")
# Remove low abundance orders
pyro = filter_taxa(pyro, function(x) sum(x) > 1, TRUE)
taxabarplot <- plot_bar(pyro, x= "SampleID", fill = "Order") + scale_y_continuous(expand = c(0, 0)) + ggtitle("") + theme(legend.title=element_blank()) + geom_bar(aes(fill=Order), stat="identity", position="stack", width =0.9) +theme_classic() + theme(text = element_text(size=14))+theme(axis.text.x = element_text(angle = 90)) + xlab("Sample") + ylab("Relative Abundance (%)") + theme(text = element_text(size=14)) + scale_fill_manual(values= pp_palette)+ xlab("")
taxabarplotBar plot at family level:
ps_nochloro_RA<- transform_sample_counts(ps_nochloro, function(OTU) 100* OTU/sum(OTU))
#Group OTU by family and sum abundances
ps_nochloro_RA_glomF<- tax_glom(ps_nochloro_RA, 'Family')
#Keep only the Pyro samples
pyro <- subset_samples(ps_nochloro_RA_glomF, Project == "Pyro")
#Remove low abundance orders
pyro = filter_taxa(pyro, function(x) sum(x) > 1, TRUE)
taxabarplot<-plot_bar(pyro, x= "SampleID", fill = "Family") + scale_y_continuous(expand = c(0, 0)) + ggtitle("") + theme(legend.title=element_blank()) + geom_bar(aes(fill=Family), stat="identity", position="stack", width =0.9) +theme_classic() + theme(text = element_text(size=14))+theme(axis.text.x = element_text(angle = 90)) + xlab("Sample") + ylab("Relative Abundance (%)") + theme(text = element_text(size=14)) + scale_fill_manual(values= pp_palette)+ xlab("")
taxabarplotParse data:
#Filter SampleID
pyroE <- subset_samples(pyro, !(SampleID %in% c("L14_6283", "L15_6283", "L17_6283", "L18_6283")))
#Remove low abundance orders
pyroE = filter_taxa(pyroE, function(x) sum(x) > 1, TRUE)
taxabarplot<-plot_bar(pyroE, x= "SampleID", fill = "Family") + scale_y_continuous(expand = c(0, 0)) + ggtitle("") + theme(legend.title=element_blank()) + geom_bar(aes(fill=Family), stat="identity", position="stack", width =0.9) +theme_classic() + theme(text = element_text(size=14))+theme(axis.text.x = element_text(angle = 90)) + xlab("Sample") + ylab("Relative Abundance (%)") + theme(text = element_text(size=14)) + scale_fill_manual(values= pp_palette)+ xlab("")
taxabarplotDetermine Unique OTU’s and abundances in P. bahamense data:
## [1] "12a986d7ed4f3be699daa0feaba41d8a" "07ae9a3ebf44d4a821d001debe7d9b0c"
## [3] "e9034a72595933ebd5c5df08eb4a08fc"
12a986d7ed4f3be699daa0feaba41d8a GCAGCAGTGGGGAATATTGGACAATGGGCGAAAGCCTGATCCAGCCATGCCGCGTGAGTGACGAAGGCCTTAGGGTTGTAAAACTCTTTTGGCGGGGACGATAATGACGGTACCCGCAGAATAAGCCCCGGCTAACTTCGTGCCAGCAGCCGCGGTAATACGAAGGGGGCTAGCGTTGTTCGGAATCACTGGGCGTAAAGCGCACGTAGGCGGACTGGTCAGTTGGGGGTGAAATCCCAGGGCTCAACCCTGGAACTGCCTCCAATACTGCCAGTCTTGAGTCCGAGAGAGGTGAGTGGAATTCCTAGTGTAGAGGTGAAATTCGTAGATATTAGGAAGAACACCAGTGGCGAAGGCGGCTCACTGGCTCGGTACTGACGCTGAGGTGCGAAAGCGTGGGGAGCAAACAGGATTAGATAC
uncultured sanger sequence from Mangrove Sediment in Hong Kong - GenBank: MH091208.1, MH091199.1 (2 best hits, 100% ) Diversity and dynamics of microbial community structure in different mangrove, marine and freshwater sediments during anaerobic debromination of PBDEs
Diversity of salt marsh prokaryotes, M.A. Moran, uncultured Hyphomicrobiaceae bacterium - lon=81.2797W, lat=31.3884N; sediment 14-16cm collected on Feb 01, 2002, Sapelo Island Microbial Observatory Dean Creek Marsh sampling site” (3rd best hit 99%)
pE_m %>% filter(OTU == "12a986d7ed4f3be699daa0feaba41d8a") %>% summarize(mean = mean(Abundance), min = min(Abundance), max = max(Abundance))## mean min max
## 1 40.82156 27.43475 50.63596
e0169d54945a8a584037e96365ec7bb3 GCTGCAGTGGGGAATATTGGACAATGGGCGAAAGCCTGATCCAGCCATGCCGCGTGAGTGACGAAGGCCTTAGGGTTGTAAAACTCTTTTGGCGGGGACGATAATGACGGTACCCGCAGAATAAGCCCCGGCTAACTTCGTGCCAGCAGCCGCGGTAATACGAAGGGGGCTAGCGTTGTTCGGAATCACTGGGCGTAAAGCGCACGTAGGCGGACTGGTCAGTTGGGGGTGAAATCCCAGGGCTCAACCCTGGAACTGCCTCCAATACTGCCAGTCTTGAGTCCGAGAGAGGTGAGTGGAATTCCTAGTGTAGAGGTGAAATTCGTAGATATTAGGAAGAACACCAGTGGCGAAGGCGGCTCACTGGCTCGGTACTGACGCTGAGGTGCGAAAGCGTGGGGAGCAAACAGGATTAGATAC
MH091208.1 uncultured bacteria - Diversity and dynamics of microbial community structure in different mangrove, marine and freshwater sediments during anaerobic debromination of PBDEs - isolation_source=“mangrove sediments” in Hong Kong, 99%
pE_m %>% filter(OTU == "e0169d54945a8a584037e96365ec7bb3") %>% summarize(mean = mean(Abundance), min = min(Abundance), max = max(Abundance))## mean min max
## 1 NaN Inf -Inf
aa0c627da86b1ddbba980c9dfd6ac380 GCAGCAGTGGGGAATATTGGACAATGGGCGAAAGCCTGATCCAGCCATGCCGCGTGAGTGACGAAGGCCTTAGGGTTGTAAAACTCTTTTGGCGGGGACGATAATGACGGTACCCGCAGAATAAGCCCCGGCTAACTTCGTGCCAGCAGCCGCGGTAATACGAAGGGGGCTAGCGTTGTTCGGAATCACTGGGCGTAAAGCGCACGTAGGCGGACTGGTCAGTTGGGGGTGAAATCCCAGGGCTCAACCCTGGAACTGCCTCCAATACTGCCAGTCTTGAGTCCGAGAGAGGTGAGTGGAATTCCTAGTGTAGAGGTGAAATTCGTAGATATTAGGAAGAACACCAGTGGCGAAGGCGACTCACTGGCTCGGTACTGACGCTGAGGTGCGAAAGCGTGGGGAGCAAACAGGATTAGATAC
MH091208.1 uncultured bacteria - Diversity and dynamics of microbial community structure in different mangrove, marine and freshwater sediments during anaerobic debromination of PBDEs - Hong Kong isolation_source=“mangrove sediments” , 99%
pE_m %>% filter(OTU == "aa0c627da86b1ddbba980c9dfd6ac380") %>% summarize(mean = mean(Abundance), min = min(Abundance), max = max(Abundance))## mean min max
## 1 NaN Inf -Inf
o__Rhizobiales f__Hyphomicrobiaceae g__Filomicrobium s__uncultured_Alphaproteobacteria
f6a0f1f2d2a083164d1fb145b000b96b GCAGCAGTGGGGAATATTGGACAATGGGCGAAAGCCTGATCCAGCCATGCCGCGTGAGTGACGAAGGCCTTAGGGTTGTAAAACTCTTTTGGTGGGGACGATAATGACGGTACCCGCAGAATAAGCCCCGGCTAACTTCGTGCCAGCAGCCGCGGTAATACGAAGGGGGCTAGCGTTGTTCGGAATCACTGGGCGTAAAGCGCACGTAGGCGGACTGGTCAGTTGGGGGTGAAATCCCAGGGCTCAACCCTGGAACTGCCTCCAATACTGCCAGTCTTGAGTCCGAGAGAGGTGAGTGGAATTCCTAGTGTAGAGGTGAAATTCGTAGATATTAGGAAGAACACCAGTGGCGAAGGCGGCTCACTGGCTCGGTACTGACGCTGAGGTGCGAAAGCGTGGGGAGCAAACAGGATTAGATAC
MH091208.1 Uncultured bacteria Diversity and dynamics of microbial community structure in different mangrove, marine and freshwater sediments during anaerobic debromination of PBDEs - isolation_source=“mangrove sediments” - 99%
pE_m %>% filter(OTU == "f6a0f1f2d2a083164d1fb145b000b96b") %>% summarize(mean = mean(Abundance), min = min(Abundance), max = max(Abundance))## mean min max
## 1 NaN Inf -Inf
ff08ec2d8a67e142bc459567275bb888 GCAGCAGTGGGGAATATTGGACAATGGGCGAAAGCCTGATCCAGCCATGCCGCGTGAGTGACGAAGGCCTTAGGGTTGTAAAACTCTTTTGGCGGGGACGATAATGACGGTACCCGCAGAATAAGCCCCGGCTAACTTCGTGCCAGCAGCCGCGGTAATACGAAGGGGGCTAGCGTTGTTCGGAATCACTGGGCGTAAAGCGCACGTAGGCGGACTGGTCAGTTGGGGGTGAAATCCCAGGGCTCAACCCTGGAACTGCCTCCAATACTGCCAGTCTTGAGTCCGAGAGAGGTGACTGGAATTCCTAGTGTAGAGGTGAAATTCGTAGATATTAGGAAGAACACCAGTGGCGAAGGCGGCTCACTGGCTCGGTACTGACGCTGAGGTGCGAAAGCGTGGGGAGCAAACAGGATTAGATAC
MH091208.1 uncultured bacteria Diversity and dynamics of microbial community structure in different mangrove, marine and freshwater sediments during anaerobic debromination of PBDEs Hong Kong isolation_source=“mangrove sediments” 99%
EU488009.1 uncultured bacteria Characterization of the lucinid bivalve-bacteria symbiotic system: the significance of the geochemical habitat on bacterial symbiont diversity and phylogeny isolation_source=“siliciclastic sedment from Thalassia sea grass bed” 99%
other marine sediments, estuary sediments
pE_m %>% filter(OTU == "ff08ec2d8a67e142bc459567275bb888") %>% summarize(mean = mean(Abundance), min = min(Abundance), max = max(Abundance))## mean min max
## 1 NaN Inf -Inf
07ae9a3ebf44d4a821d001debe7d9b0c GCAGCAGTGGGGAATCTTGCGCAATGGGCGAAAGCCTGACGCAGCGACGCCGCGTGCGGGAAGACGGCCTTCGGGTTGTAAACCGCTTTCAGCAGGGACGAAATTGACGGTACCTGCAGAAGAAGCTCCGGCCAACTACGTGCCAGCAGCCGCGGTAAGACGTAGGGGGCGAGCGTTGTCCGGAATCATTGGGCGTAAAGGGCTCCTAGGTGGTTCAGTAAGTCGACTGTGAAAATCCAAGGCTCAACCTTGGGACGCCAGTCGATACTGCTGTGACTCGAGTTCGGTAGAGGAGTGTGGAATTCCTGGTGTAGCGGTGAAATGCGCAGATATCAGGAGGAACACCAACGGCGAAGGCAGCACTCTGGGCCGATACTGACACTGAAGAGCGAAAGCGTGGGGAGCAAACAGGATTAGATAC
KM840989.1 uncultured bacteria - Microbial diversity of indigenous bacteria in a 129I contaminated groundwater plume at the Hanford Site, Washington - Sanger, isolation_source=“iodine (I-129) contaminated groundwater” 100%
HF558551.1 uncultured bacteria - Iron- and Sulphur- cycling bacteria mobilize copper in a multiple extreme mine tailings in the Atacama Desert, Chile - isolation_source=“tailing material” - 100%
JQ427269.1 uncultured bacteria - Bacterial diversity in an alkaline saline soil spiked with anthracene, samger sequenced, isolation_source=“soil”, 100%
JQ665390.1 uncultured actinobacterium, Diversity of unculturable Actinomycetes in coastal wetlands of the Yellow River estuary, isolation_source=“soil sample from coastal wetlands”, 99%
pE_m %>% filter(OTU == "07ae9a3ebf44d4a821d001debe7d9b0c") %>% summarize(mean = mean(Abundance), min = min(Abundance), max = max(Abundance))## mean min max
## 1 33.68851 22.47499 44.42002
9b8accad19c30d7a69a6290553111166 GCTGCAGTGGGGAATCTTGCGCAATGGGCGAAAGCCTGACGCAGCGACGCCGCGTGCGGGAAGACGGCCTTCGGGTTGTAAACCGCTTTCAGCAGGGACGAAATTGACGGTACCTGCAGAAGAAGCTCCGGCCAACTACGTGCCAGCAGCCGCGGTAAGACGTAGGGGGCGAGCGTTGTCCGGAATCATTGGGCGTAAAGGGCTCCTAGGTGGTTCAGTAAGTCGACTGTGAAAATCCAAGGCTCAACCTTGGGACGCCAGTCGATACTGCTGTGACTCGAGTTCGGTAGAGGAGTGTGGAATTCCTGGTGTAGCGGTGAAATGCGCAGATATCAGGAGGAACACCAACGGCGAAGGCAGCACTCTGGGCCGATACTGACACTGAAGAGCGAAAGCGTGGGGAGCAAACAGGATTAGATAC
KM840989.1 uncultured bacteria - Microbial diversity of indigenous bacteria in a 129I contaminated groundwater plume at the Hanford Site, Washington - isolation_source=“iodine (I-129) contaminated groundwater” - 99%
pE_m %>% filter(OTU == "9b8accad19c30d7a69a6290553111166") %>% summarize(mean = mean(Abundance), min = min(Abundance), max = max(Abundance))## mean min max
## 1 NaN Inf -Inf
e9034a72595933ebd5c5df08eb4a08fc GCAGCAGTGGGGAATCTTGGACAATGGGGGCAACCCTGATCCAGCCATGCCGCGTGTGTGAAGAAGGCCTTCGGGTTGTAAAGCACTTTCAGCAGGGAGGAAGGCTTACCCCTAATACGGGTGAGTACTTGACGTTACCTGCAGAAGAAGCACCGGCTAATTTCGTGCCAGCAGCCGCGGTAATACGAAAGGTGCAAGCGTTAATCGGAATTACTGGGCGTAAAGCGCGCGTAGGCGGTGTGTTAAGTCGGATGTGAAAGCCCAGGGCTCAACCTTGGAATTGCATCCGATACTGGCACGCTAGAGTGCAGTAGAGGGAGGTGGAATTTCCGGTGTAGCGGTGAAATGCGTAGAGATCGGAAGGAACACCAGTGGCGAAGGCGGCCTCCTGGACTGACACTGACGCTGAGGTGCGAAAGCGTGGGGAGCAAACAGGATTAGATAC
PP516259.1 organism=“Alloalcanivorax venustensis” Exploration and conservation of bacterial community from the Arabian Sea seamount isolation_source=“Arabian Sea water”, 100%
pE_m %>% filter(OTU == "e9034a72595933ebd5c5df08eb4a08fc") %>% summarize(mean = mean(Abundance), min = min(Abundance), max = max(Abundance))## mean min max
## 1 25.4058 18.12139 34.43303
c02844d75d512c1510de7334b6687731 GCAGCAGTGGGGAATATTGGACAATGGGCGAAAGCCTGATCCAGCCATGCCGCGTGAGTGACGAAGGCCTTAGGGTTGTAAAGCTCTTTTGGCGGGGAAGATAATGACGGTACCCGCAGAATAAGCTCCGGCTAACTTCGTGCCAGCAGCCGCGGTAATACGAAGGGAGCTAGCGTTGTTCGGAATCACTGGGCGTAAAGCGCACGTAGGCGGATTTGTTAGTCAGGGGTGAAATCCCGGGGCTCAACCCCGGAACTGCCTTTGATACTGCAAATCTCGAGTCCGAGAGAGGTGGGTGGAATTCCTAGTGTAGAGGTGAAATTCGTAGATATTAGGAAGAACACCGGTGGCGAAGGCGGCCCACTGGCTCGGTACTGACGCTGAGGTGCGAAAGCGTGGGGAGCAAACAGGATTAGATAC
OK235762.1 Uncultured Pedomicrobium sp. isolation_source=“date palm rhizosphere soil” Saudi Arabia - 99%
HQ697761.1 Uncultured bacteria isolation_source=“hydrocarbon contaminated saline-alkali soil” China - 99%
99% to several other contaminated soil sites (real and experimental)
KP098952.1 Uncultured bacteria - The microbiome of methanol-utilizing denitrification systems contains new bacterial groups - isolation_source=“denitrification bioreactor” - 99%
pE_m %>% filter(OTU == "c02844d75d512c1510de7334b6687731") %>% summarize(mean = mean(Abundance), min = min(Abundance), max = max(Abundance))## mean min max
## 1 NaN Inf -Inf
Plot shannon diversity and relative abundance:
# -------------------------------
# Subset Pyro samples and filter
# -------------------------------
pyro <- subset_samples(ps_nochloro_RA_glomO, Project == "Pyro")
pyro <- subset_samples(pyro, Treatment != "AB")
pyro <- filter_taxa(pyro, function(x) sum(x) > 1, TRUE)
pyro <- subset_samples(pyro,
Treatment %in% c("1003", "Replete") &
TreatRep %in% c("1003_1", "1_2"))
# Set factor levels for Notes_Filter
pyro@sam_data$Notes_Filter <- factor(as.character(pyro@sam_data$Notes_Filter),
levels = c("10", "2"))
# Clean Order names in tax_table
tax_table(pyro)[, "Order"] <- gsub(" o__", "", tax_table(pyro)[, "Order"])
# -------------------------------
# Plot relative abundance barplot
# -------------------------------
barplot <- plot_bar(pyro, x = "TreatRep", fill = "Order") +
geom_bar(aes(fill = Order), stat = "identity", position = "stack", width = 1) +
scale_y_continuous(expand = c(0, 0)) +
scale_fill_manual(values = pp_palette) +
facet_grid(Notes_Filter ~ Treatment, scales = "free") +
scale_x_discrete(
expand = c(0.5, 0.05),
labels = c(
"1003_1" = "IRL",
"1_2" = "OTB"
)
) +
theme_classic() +
theme(
legend.title = element_blank(),
legend.position = "right",
legend.justification = "center",
legend.box = "horizontal",
axis.text.x = element_text(angle = 0, vjust = 0.5, size = 16, color = "black"),
axis.text.y = element_text(size = 20, color = "black"),
axis.title.x = element_text(size = 20, color = "black"),
axis.title.y = element_text(size = 20, color = "black"),
axis.ticks.x = element_line(color = "black", linewidth = 0.9),
axis.ticks.y = element_line(color = "black", linewidth = 0.9),
strip.text.x = element_blank(), # remove top facet labels
strip.text.y = element_text(size = 20, color = "black"), # row facets
text = element_text(size = 20, color = "black")
) +
ylab(expression("Relative Abundance (%)")) +
xlab("")
barplot# -------------------------------
# Calculate Shannon index
# -------------------------------
pyro <- subset_samples(ps_nochloro_RA_glomO, Project == "Pyro")
pyro <- subset_samples(pyro, TreatRep %in% c("1003_1", "1_2", "1_1", "1_3"))
pyro <- subset_samples(pyro, Treatment != "AB")
shannon_df <- estimate_richness(pyro, measures = "Shannon")
# Extract metadata and combine with diversity
sample_meta <- as(sample_data(pyro), "data.frame")
shannon_df <- cbind(shannon_df, sample_meta)
# -------------------------------
# Plot Shannon diversity
# -------------------------------
shannon_plot <- ggplot(shannon_df, aes(x = Treatment, y = Shannon, color = Notes_Filter)) +
# Color mapping
scale_color_manual(values = c("10" = "navy", "2" = "maroon")) +
# Points
geom_point(size = 5) +
# Custom x-axis labels
scale_x_discrete(labels = c(
"1003" = "IRL",
"Replete" = "OTB"
)) +
# Classic theme
theme_classic() +
# Axis labels
ylab("Shannon Diversity Index") +
xlab("") +
labs(color = "Filter Size") +
# Theme adjustments
theme(
legend.title = element_text(size = 20),
axis.text.x = element_text(angle = 0, vjust = 0.5, size = 16, color = "black"),
axis.text.y = element_text(size = 20, color = "black"),
axis.title.x = element_text(size = 20, color = "black"),
axis.title.y = element_text(size = 20, color = "black"),
axis.ticks.x = element_line(color = "black", linewidth = 0.9),
axis.ticks.y = element_line(color = "black", linewidth = 0.9),
strip.text.x = element_blank(), # remove top facet labels
strip.text.y = element_text(size = 20, color = "black"), # row facets
text = element_text(size = 20, color = "black"))
# Print Shannon plot
shannon_plotChlorophyll-a fluorescence (RFU) data from P. bahamense cultures were used to visualize growth patterns, quantify maximum biomass, and estimate specific growth rates under media conditions deficient in media lacking one targeted B vitamin. Specific growth rates were computed using the “alleasy_linear” function of the “growthrates” package in R Studio (Hall et al., 2014; RStudio Team, 2020).
Import the Dataset:
fluor<-as.data.frame(read_csv("fluor.csv"))
fluor$Time <- as.numeric(sub('.', '', fluor$Time))
fluor$Date <- mdy(fluor$Date)
fluor$Time <- fluor$Time * 48 /24LD_OTB_all_growth <- ggplot(
fluor %>%
filter(
Treatment %in% c("-Cobalamin", "Replete", "-Thiamine", "-Biotin"),
Round %in% c(1, 2,3,4)),
aes(x = Time, y = `Chlorophyll-A (RFU)`, color = Treatment)) +
geom_point(size = 3) + geom_smooth(se=FALSE)+
# Manual color scale
scale_color_manual(
values = c(
"-Cobalamin" = "#102E50",
"Replete" = "#F5C45E",
"-Thiamine" = "deeppink4",
"-Biotin" = "#667028"
),
labels = c(
"-Cobalamin" = expression(" —B"[12]),
"-Thiamine" = expression("—thiamine"),
"-Biotin" = expression("—biotin"),
"Replete" = "replete"),
name = "Treatment") +
labs(
x = "Time (days)",
y = expression(paste("Chlorophyll-", italic("a"), " fluorescence (RFU)"))
) +
facet_grid(~Round) +
theme_minimal() +
theme(
panel.border = element_rect(fill = NA, color = "black"),
panel.grid = element_blank(),
axis.title = element_text(size = 20),
axis.text = element_text(size = 20, color = "black"),
strip.text = element_text(size = 20, color = "black"),
legend.title = element_text(size = 18, color = "black"),
legend.text = element_text(size = 18, color = "black"),
legend.position = "bottom",
axis.ticks = element_line(color = "black", linewidth = 0.5),
axis.ticks.length = unit(4, "pt")
)
LD_OTB_all_growthmax_bio <- fluor %>%
filter(Round %in% c(1, 2, 3, 4)) %>%
filter(Treatment %in% c("-Cobalamin", "Replete", "-Biotin", "-Thiamine")) %>%
group_by(Treatment, Round, Rep) %>%
summarise(Max_Chl = max(`Chlorophyll-A (RFU)`, na.rm = TRUE), .groups = "drop")
print(max_bio)## # A tibble: 48 × 4
## Treatment Round Rep Max_Chl
## <chr> <dbl> <dbl> <dbl>
## 1 -Biotin 1 1 2911.
## 2 -Biotin 1 2 3220.
## 3 -Biotin 1 3 6239.
## 4 -Biotin 2 1 6771.
## 5 -Biotin 2 2 4230.
## 6 -Biotin 2 3 3492.
## 7 -Biotin 3 1 4288.
## 8 -Biotin 3 2 3237.
## 9 -Biotin 3 3 1767.
## 10 -Biotin 4 1 2049.
## # ℹ 38 more rows
Run an ANOVA on maximum RFU values and Tukey HSD post-hoc test for treatment
## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = Max_Chl ~ Treatment, data = max_bio)
##
## $Treatment
## diff lwr upr p adj
## -Cobalamin--Biotin -3365.871 -4602.055 -2129.687 0.0000000
## -Thiamine--Biotin 1961.220 725.036 3197.404 0.0006407
## Replete--Biotin -757.515 -1993.699 478.669 0.3694455
## -Thiamine--Cobalamin 5327.091 4090.907 6563.275 0.0000000
## Replete--Cobalamin 2608.356 1372.172 3844.540 0.0000068
## Replete--Thiamine -2718.735 -3954.919 -1482.551 0.0000030
## Df Sum Sq Mean Sq F value Pr(>F)
## Treatment 3 174966709 58322236 45.35 1.63e-13 ***
## Residuals 44 56590796 1286154
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Prepare data for plotting maximum RFU values:
# Calculate mean and SD
summary_stats <- max_bio %>%
group_by(Treatment, Round) %>%
summarise(
mean_chl = mean(Max_Chl, na.rm = TRUE),
sd_chl = sd(Max_Chl, na.rm = TRUE),
.groups = "drop"
)
# Make Treatment a factor in both datasets
max_bio$Treatment <- factor(max_bio$Treatment,
levels = c( "-Biotin", "-Cobalamin","-Thiamine", "Replete"))
summary_stats$Treatment <- factor(summary_stats$Treatment,
levels = c("-Cobalamin", "-Thiamine", "-Biotin", "Replete"))
max_bio$Treatment <- factor(
max_bio$Treatment,
levels = c("-Cobalamin", "Replete", "-Thiamine", "-Biotin")
)
summary_stats$Treatment <- factor(
summary_stats$Treatment,
levels = c("-Cobalamin", "Replete", "-Thiamine", "-Biotin"))Plot maximum RFU values:
max_bio_LD <- ggplot(max_bio, aes(x = Treatment, y = Max_Chl, color = Treatment)) +
# Raw points (per round)
geom_point(size = 3, shape = 16, position = position_jitter(width = 0.1)) +
# Error bars: mean ± SD across replicates
geom_errorbar(
data = summary_stats,
aes(
x = Treatment,
y = mean_chl,
ymin = mean_chl - sd_chl,
ymax = mean_chl + sd_chl,
color = Treatment
),
width = 0.3,
position = position_dodge(width = 0.6),
inherit.aes = FALSE
) +
# Open circles for replicate-averaged means
geom_point(
data = summary_stats,
aes(x = Treatment, y = mean_chl, color = Treatment),
shape = 21,
size = 5,
stroke = 1.2,
position = position_dodge(width = 0.6),
inherit.aes = FALSE
) +
# Labels
labs(
x = "Media Treatment",
y = expression(paste("Maximum RFU"))
) +
# Manual colors
scale_color_manual(
values = c(
"-Cobalamin" = "#102E50",
"-Thiamine" = "deeppink4",
"-Biotin" = "#667028",
"Replete" = "#F5C45E"
),
labels = c(
"-Cobalamin" = expression("—B"[12]),
"-Thiamine" = expression("—Thiamine"),
"-Biotin" = expression("—Biotin"),
"Replete" = "Replete"
),
name = "Treatment",
guide = guide_legend(
override.aes = list(
shape = 21,
fill = "white",
stroke = 1.2,
size = 2
)
)
) +
# X-axis label formatting
scale_x_discrete(
labels = c(
"-Biotin" = expression("\u2013Biotin"),
"-Cobalamin" = expression("\u2013B"[12]),
"-Thiamine" = expression("\u2013Thiamine"),
"Replete" = "Replete"
)
) +
scale_y_continuous(labels = comma) +
# Theme
theme_minimal() +
theme(
panel.border = element_rect(fill = NA, color = "black"),
panel.grid = element_blank(),
axis.text.y = element_text(size = 20, color = "black"),
axis.text.x = element_text(angle = 45, hjust = 1, size = 16, color = "black"),
axis.title.y = element_text(size = 22),
axis.title.x = element_text(size = 20, margin = margin(t = 10)),
strip.text.x = element_text(size = 18),
legend.title = element_text( size = 20, color = "black"),
legend.text = element_text( size = 16, color = "black"),
legend.position = "none",
axis.ticks = element_line(color = "black", linewidth = 0.5),
axis.ticks.length = unit(4, "pt")
) +
guides(color = "none", fill = "none", shape = "none")+
coord_cartesian(ylim = c(0, 10000))
# Print plot
max_bio_LDSpecific Growth Rates Biotin trial 1 using all_easylinear:
bio1 <- fluor %>%
filter(Treatment == "-Biotin", Round %in% c(1) )%>%
rename(Chl = `Chlorophyll-A (RFU)`) %>%
droplevels()
bio_lin1 <- all_easylinear(Chl ~ Time|Round+Treatment+Rep, bio1, h = 7, quota = 1)
par(mfrow = c(3, 3))
par(mar = c(2, 2, 2, 1))
plot(bio_lin1, log = "y")
summary(bio_lin1)## $`1:-Biotin:1`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## -0.1426 -0.3358 0.3733 0.3348 0.1439 -0.2512 -0.1225
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.66156 0.21185 17.284 1.19e-05 ***
## x 0.20749 0.02938 7.063 0.00088 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3109 on 5 degrees of freedom
## Multiple R-squared: 0.9089, Adjusted R-squared: 0.8907
## F-statistic: 49.88 on 1 and 5 DF, p-value: 0.0008799
##
##
## $`1:-Biotin:2`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## -0.2611 0.1266 0.1952 0.1704 -0.1053 -0.1483 0.0224
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.77140 0.13104 28.78 9.49e-07 ***
## x 0.21617 0.01817 11.90 7.40e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1923 on 5 degrees of freedom
## Multiple R-squared: 0.9659, Adjusted R-squared: 0.959
## F-statistic: 141.5 on 1 and 5 DF, p-value: 7.396e-05
##
##
## $`1:-Biotin:3`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## -0.018655 -0.066394 0.061972 0.130916 -0.079399 -0.037937 0.009497
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.479439 0.142763 17.37 1.16e-05 ***
## x 0.257843 0.007742 33.30 4.59e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.08194 on 5 degrees of freedom
## Multiple R-squared: 0.9955, Adjusted R-squared: 0.9946
## F-statistic: 1109 on 1 and 5 DF, p-value: 4.589e-07
## y0 y0_lm mumax lag
## 1:-Biotin:1 33.75 38.92188 0.2074863 -0.6871588
## 1:-Biotin:2 33.46 43.44102 0.2161749 -1.2076036
## 1:-Biotin:3 49.69 11.93456 0.2578427 5.5319204
## 1:-Biotin:1.r2 1:-Biotin:2.r2 1:-Biotin:3.r2
## 0.9088927 0.9658730 0.9955119
Specific Growth Rates Biotin trial 2 using all_easylinear
bio2 <- fluor %>%
filter(Treatment == "-Biotin", Round %in% c(2) )%>%
rename(Chl = `Chlorophyll-A (RFU)`) %>%
droplevels()
bio_lin2 <- all_easylinear(Chl ~ Time|Round+Treatment+Rep, bio2, h = 7, quota = 1)
par(mfrow = c(3, 3))
par(mar = c(2, 2, 2, 1))
plot(bio_lin2, log = "y")
summary(bio_lin2)## $`2:-Biotin:1`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## -1.29805 0.90387 0.62634 0.53091 -0.40115 -0.02684 -0.33508
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.63417 1.14283 1.430 0.21213
## x 0.37606 0.07849 4.791 0.00492 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.8307 on 5 degrees of freedom
## Multiple R-squared: 0.8211, Adjusted R-squared: 0.7854
## F-statistic: 22.96 on 1 and 5 DF, p-value: 0.004922
##
##
## $`2:-Biotin:2`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## -0.06137 -0.08479 -0.10821 0.37956 -0.09211 0.27054 -0.30362
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.36851 0.26561 12.682 5.42e-05 ***
## x 0.24466 0.02466 9.921 0.000178 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.261 on 5 degrees of freedom
## Multiple R-squared: 0.9517, Adjusted R-squared: 0.942
## F-statistic: 98.42 on 1 and 5 DF, p-value: 0.0001776
##
##
## $`2:-Biotin:3`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## 0.24068 -0.16024 -0.25709 -0.01833 0.16293 0.11462 -0.08257
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.10851 0.16756 24.52 2.10e-06 ***
## x 0.29515 0.01873 15.76 1.87e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1983 on 5 degrees of freedom
## Multiple R-squared: 0.9803, Adjusted R-squared: 0.9763
## F-statistic: 248.2 on 1 and 5 DF, p-value: 1.873e-05
## y0 y0_lm mumax lag
## 2:-Biotin:1 98.46 5.125203 0.3760631 7.859002
## 2:-Biotin:2 96.62 29.035114 0.2446627 4.914031
## 2:-Biotin:3 107.23 60.856127 0.2951547 1.919209
## 2:-Biotin:1.r2 2:-Biotin:2.r2 2:-Biotin:3.r2
## 0.8211458 0.9516543 0.9802552
Specific Growth Rates Biotin trial 3 using all_easylinear:
bio3 <- fluor %>%
filter(Treatment == "-Biotin", Round %in% c(3) )%>%
rename(Chl = `Chlorophyll-A (RFU)`) %>%
droplevels()
bio_lin3 <- all_easylinear(Chl ~ Time|Round+Treatment+Rep, bio3, h = 7, quota = 1)
par(mfrow = c(3, 3))
par(mar = c(2, 2, 2, 1))
plot(bio_lin3, log = "y")
summary(bio_lin3)## $`3:-Biotin:1`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## -0.07856 -0.22592 0.15259 0.26659 0.21258 -0.23432 -0.09296
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.81040 0.35557 7.904 0.000522 ***
## x 0.23963 0.02156 11.115 0.000103 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2282 on 5 degrees of freedom
## Multiple R-squared: 0.9611, Adjusted R-squared: 0.9533
## F-statistic: 123.5 on 1 and 5 DF, p-value: 0.0001028
##
##
## $`3:-Biotin:2`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## 0.2353339 -0.1112728 -0.0556536 -0.2798767 0.0004855 0.2056337 0.0053499
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.96414 0.23391 12.67 5.44e-05 ***
## x 0.26925 0.01849 14.56 2.76e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1957 on 5 degrees of freedom
## Multiple R-squared: 0.977, Adjusted R-squared: 0.9723
## F-statistic: 212 on 1 and 5 DF, p-value: 2.76e-05
##
##
## $`3:-Biotin:3`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## 0.01417 -0.22670 0.78580 0.49481 -1.90075 0.22234 0.61033
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.56645 2.28724 0.248 0.814
## x 0.18591 0.09401 1.978 0.105
##
## Residual standard error: 0.9949 on 5 degrees of freedom
## Multiple R-squared: 0.4389, Adjusted R-squared: 0.3267
## F-statistic: 3.911 on 1 and 5 DF, p-value: 0.1049
## y0 y0_lm mumax lag
## 3:-Biotin:1 98.82 16.616595 0.2396343 7.440080
## 3:-Biotin:2 63.40 19.378045 0.2692450 4.402396
## 3:-Biotin:3 82.43 1.761997 0.1859051 20.685290
## 3:-Biotin:1.r2 3:-Biotin:2.r2 3:-Biotin:3.r2
## 0.9611016 0.9769567 0.4388915
Specific Growth Rates Biotin trial 4 using all_easylinear
bio4 <- fluor %>%
filter(Treatment == "-Biotin", Round %in% c(4) )%>%
rename(Chl = `Chlorophyll-A (RFU)`) %>%
droplevels()
bio_lin4 <- all_easylinear(Chl ~ Time|Round+Treatment+Rep, bio4, h = 7, quota = 1)
par(mfrow = c(3, 3))
par(mar = c(2, 2, 2, 1))
plot(bio_lin4, log = "y")
summary(bio_lin4)## $`4:-Biotin:1`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## 0.30806 0.04318 -0.14397 -0.50617 -0.10602 0.24216 0.16276
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.26650 0.53423 2.371 0.063898 .
## x 0.23271 0.02897 8.032 0.000484 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3066 on 5 degrees of freedom
## Multiple R-squared: 0.9281, Adjusted R-squared: 0.9137
## F-statistic: 64.51 on 1 and 5 DF, p-value: 0.0004838
##
##
## $`4:-Biotin:2`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## -0.19446 0.54358 -0.06320 -0.20050 -0.37665 0.05642 0.23480
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.80368 0.59222 3.046 0.02857 *
## x 0.21685 0.03212 6.752 0.00108 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3399 on 5 degrees of freedom
## Multiple R-squared: 0.9012, Adjusted R-squared: 0.8814
## F-statistic: 45.58 on 1 and 5 DF, p-value: 0.001082
##
##
## $`4:-Biotin:3`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## 0.02364 0.05129 0.03152 -0.14500 -0.02402 -0.04135 0.10391
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.488502 0.120401 20.67 4.91e-06 ***
## x 0.233579 0.008269 28.25 1.04e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.08751 on 5 degrees of freedom
## Multiple R-squared: 0.9938, Adjusted R-squared: 0.9925
## F-statistic: 797.9 on 1 and 5 DF, p-value: 1.041e-06
## y0 y0_lm mumax lag
## 4:-Biotin:1 57.52 3.548409 0.2327066 11.970583
## 4:-Biotin:2 67.67 6.071954 0.2168469 11.118270
## 4:-Biotin:3 97.55 12.043218 0.2335787 8.955711
## 4:-Biotin:1.r2 4:-Biotin:2.r2 4:-Biotin:3.r2
## 0.9280693 0.9011555 0.9937725
Specific Growth Rates Thiamine trial 1 using all_easylinear:
thia1 <- fluor %>%
filter(Treatment == "-Thiamine", Round %in% c(1) )%>%
rename(Chl = `Chlorophyll-A (RFU)`) %>%
droplevels()
thia_lin1 <- all_easylinear(Chl ~ Time|Round+Treatment+Rep, thia1, h = 7, quota = 1)
par(mfrow = c(3, 3))
par(mar = c(2, 2, 2, 1))
plot(thia_lin1, log = "y")
summary(thia_lin1)## $`1:-Thiamine:1`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## 0.30173 -0.51432 -0.04376 0.21123 0.09790 0.10677 -0.15955
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.4499 0.4668 7.391 0.000713 ***
## x 0.2417 0.0283 8.539 0.000363 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2995 on 5 degrees of freedom
## Multiple R-squared: 0.9358, Adjusted R-squared: 0.923
## F-statistic: 72.91 on 1 and 5 DF, p-value: 0.0003627
##
##
## $`1:-Thiamine:2`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## -0.31890 0.23129 0.13999 -0.02326 0.25846 -0.25014 -0.03743
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.14430 0.25181 12.49 5.84e-05 ***
## x 0.26660 0.02338 11.40 9.08e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2474 on 5 degrees of freedom
## Multiple R-squared: 0.963, Adjusted R-squared: 0.9556
## F-statistic: 130 on 1 and 5 DF, p-value: 9.079e-05
##
##
## $`1:-Thiamine:3`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## -0.25586 0.36063 0.14809 -0.29385 0.02147 -0.02179 0.04130
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.33255 0.25195 13.23 4.41e-05 ***
## x 0.26006 0.02339 11.12 0.000103 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2476 on 5 degrees of freedom
## Multiple R-squared: 0.9611, Adjusted R-squared: 0.9533
## F-statistic: 123.6 on 1 and 5 DF, p-value: 0.0001027
## y0 y0_lm mumax lag
## 1:-Thiamine:1 173.20 31.49697 0.2416543 7.0536952
## 1:-Thiamine:2 44.15 23.20350 0.2666038 2.4129054
## 1:-Thiamine:3 33.19 28.00966 0.2600576 0.6525443
## 1:-Thiamine:1.r2 1:-Thiamine:2.r2 1:-Thiamine:3.r2
## 0.9358222 0.9629722 0.9611167
Specific Growth Rates Thiamine trial 2 using all_easylinear
thia2 <- fluor %>%
filter(Treatment == "-Thiamine", Round %in% c(2) )%>%
rename(Chl = `Chlorophyll-A (RFU)`) %>%
droplevels()
thia_lin2 <- all_easylinear(Chl ~ Time|Round+Treatment+Rep, thia2, h = 9, quota = 1)
par(mfrow = c(3, 3))
par(mar = c(2, 2, 2, 1))
plot(thia_lin2, log = "y")
summary(thia_lin2)## $`2:-Thiamine:1`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.4806 -0.2537 -0.1203 0.2472 0.5625
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.63297 0.26430 17.53 4.84e-07 ***
## x 0.24127 0.02348 10.27 1.79e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3638 on 7 degrees of freedom
## Multiple R-squared: 0.9378, Adjusted R-squared: 0.9289
## F-statistic: 105.5 on 1 and 7 DF, p-value: 1.79e-05
##
##
## $`2:-Thiamine:2`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.35944 -0.06839 0.00212 0.08203 0.24602
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.80494 0.11101 43.28 9.17e-10 ***
## x 0.22326 0.01166 19.15 2.64e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1806 on 7 degrees of freedom
## Multiple R-squared: 0.9813, Adjusted R-squared: 0.9786
## F-statistic: 366.7 on 1 and 7 DF, p-value: 2.636e-07
##
##
## $`2:-Thiamine:3`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.42503 -0.12665 -0.00432 0.16254 0.40225
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.54004 0.16251 27.94 1.93e-08 ***
## x 0.25614 0.01707 15.01 1.40e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2644 on 7 degrees of freedom
## Multiple R-squared: 0.9699, Adjusted R-squared: 0.9656
## F-statistic: 225.2 on 1 and 7 DF, p-value: 1.4e-06
## y0 y0_lm mumax lag
## 2:-Thiamine:1 156.08 102.81917 0.2412673 1.7300188
## 2:-Thiamine:2 131.83 122.11205 0.2232630 0.3429772
## 2:-Thiamine:3 140.09 93.69448 0.2561405 1.5704110
## 2:-Thiamine:1.r2 2:-Thiamine:2.r2 2:-Thiamine:3.r2
## 0.9378050 0.9812701 0.9698585
Specific Growth Rates Thiamine trial 3 using all_easylinear:
thia3 <- fluor %>%
filter(Treatment == "-Thiamine", Round %in% c(3) )%>%
rename(Chl = `Chlorophyll-A (RFU)`) %>%
droplevels()
thia_lin3 <- all_easylinear(Chl ~ Time|Round+Treatment+Rep, thia3, h = 9, quota = 1)
par(mfrow = c(3, 3))
par(mar = c(2, 2, 2, 1))
plot(thia_lin3, log = "y")
summary(thia_lin3)## $`3:-Thiamine:1`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.35661 -0.21370 -0.07224 0.18665 0.48952
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.79058 0.36271 4.937 0.00168 **
## x 0.25375 0.01937 13.101 3.52e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3001 on 7 degrees of freedom
## Multiple R-squared: 0.9608, Adjusted R-squared: 0.9552
## F-statistic: 171.6 on 1 and 7 DF, p-value: 3.521e-06
##
##
## $`3:-Thiamine:2`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.49865 0.01803 0.07384 0.10434 0.40557
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.81256 0.35153 5.156 0.00132 **
## x 0.31014 0.02091 14.833 1.52e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3239 on 7 degrees of freedom
## Multiple R-squared: 0.9692, Adjusted R-squared: 0.9648
## F-statistic: 220 on 1 and 7 DF, p-value: 1.517e-06
##
##
## $`3:-Thiamine:3`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.4806 -0.1946 -0.0558 0.3160 0.5117
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.99093 0.40297 4.941 0.00167 **
## x 0.25389 0.02397 10.593 1.46e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3713 on 7 degrees of freedom
## Multiple R-squared: 0.9413, Adjusted R-squared: 0.9329
## F-statistic: 112.2 on 1 and 7 DF, p-value: 1.462e-05
## y0 y0_lm mumax lag
## 3:-Thiamine:1 49.92 5.992945 0.2537505 8.354027
## 3:-Thiamine:2 72.28 6.126119 0.3101377 7.957711
## 3:-Thiamine:3 57.00 7.322315 0.2538908 8.082708
## 3:-Thiamine:1.r2 3:-Thiamine:2.r2 3:-Thiamine:3.r2
## 0.9608132 0.9691648 0.9412777
Specific Growth Rates Thiamine trial 4 using all_easylinear
thia4 <- fluor %>%
filter(Treatment == "-Thiamine", Round %in% c(4) )%>%
rename(Chl = `Chlorophyll-A (RFU)`) %>%
droplevels()
thia_lin4 <- all_easylinear(Chl ~ Time|Round+Treatment+Rep, thia4, h = 9, quota = 1)
par(mfrow = c(3, 3))
par(mar = c(2, 2, 2, 1))
plot(thia_lin4, log = "y")
summary(thia_lin4)## $`4:-Thiamine:1`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.34954 -0.14112 0.08196 0.17775 0.22950
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.4689 0.1816 19.11 2.68e-07 ***
## x 0.2533 0.0139 18.22 3.71e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2153 on 7 degrees of freedom
## Multiple R-squared: 0.9794, Adjusted R-squared: 0.9764
## F-statistic: 332 on 1 and 7 DF, p-value: 3.71e-07
##
##
## $`4:-Thiamine:2`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.269152 -0.100510 0.007589 0.106550 0.281775
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.87449 0.21124 13.61 2.72e-06 ***
## x 0.24309 0.01256 19.35 2.46e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1946 on 7 degrees of freedom
## Multiple R-squared: 0.9816, Adjusted R-squared: 0.979
## F-statistic: 374.3 on 1 and 7 DF, p-value: 2.457e-07
##
##
## $`4:-Thiamine:3`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.7143 0.0431 0.1993 0.3349 0.4359
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.65854 0.76560 3.472 0.0104 *
## x 0.23859 0.04554 5.239 0.0012 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.7055 on 7 degrees of freedom
## Multiple R-squared: 0.7968, Adjusted R-squared: 0.7678
## F-statistic: 27.45 on 1 and 7 DF, p-value: 0.0012
## y0 y0_lm mumax lag
## 4:-Thiamine:1 161.92 32.10143 0.2532588 6.389518
## 4:-Thiamine:2 76.34 17.71644 0.2430897 6.008910
## 4:-Thiamine:3 109.65 14.27538 0.2385910 8.544988
## 4:-Thiamine:1.r2 4:-Thiamine:2.r2 4:-Thiamine:3.r2
## 0.9793538 0.9816432 0.7968191
Specific Growth Rates Replete trial 1 using all_easylinear:
rep1 <- fluor %>%
filter(Treatment == "Replete", Round %in% c(1)) %>%
rename(Chl = `Chlorophyll-A (RFU)`) %>%
droplevels()
rep_lin1 <- all_easylinear(Chl ~ Time|Round+Treatment+Rep, rep1, h = 7, quota = 1)
par(mfrow = c(3, 3))
par(mar = c(2, 2, 2, 1))
plot(rep_lin1, log = "y")
summary(rep_lin1)## $`1:Replete:1`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## -0.23658 0.13397 0.19451 -0.01732 -0.04291 0.10941 -0.14107
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.29023 0.23491 14.01 3.34e-05 ***
## x 0.21573 0.01613 13.37 4.19e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1707 on 5 degrees of freedom
## Multiple R-squared: 0.9728, Adjusted R-squared: 0.9674
## F-statistic: 178.8 on 1 and 5 DF, p-value: 4.186e-05
##
##
## $`1:Replete:2`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## -0.01471 0.09138 -0.04521 -0.16653 0.07629 0.15920 -0.10042
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.03875 0.17358 17.51 1.12e-05 ***
## x 0.22585 0.01192 18.95 7.55e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1262 on 5 degrees of freedom
## Multiple R-squared: 0.9863, Adjusted R-squared: 0.9835
## F-statistic: 358.9 on 1 and 5 DF, p-value: 7.551e-06
##
##
## $`1:Replete:3`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## 0.025927 0.009719 0.219121 -0.438641 0.040203 0.154873 -0.011203
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.29711 0.31788 10.37 0.000143 ***
## x 0.22019 0.02183 10.09 0.000164 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.231 on 5 degrees of freedom
## Multiple R-squared: 0.9531, Adjusted R-squared: 0.9438
## F-statistic: 101.7 on 1 and 5 DF, p-value: 0.0001641
## y0 y0_lm mumax lag
## 1:Replete:1 66.01 26.84907 0.2157266 4.169977
## 1:Replete:2 87.82 20.87905 0.2258538 6.360502
## 1:Replete:3 94.03 27.03446 0.2201909 5.661004
## 1:Replete:1.r2 1:Replete:2.r2 1:Replete:3.r2
## 0.9727947 0.9862601 0.9531494
Specific Growth Rates Replete trial 2 using all_easylinear:
rep2 <- fluor %>%
filter(Treatment == "Replete", Round %in% c(2)) %>%
rename(Chl = `Chlorophyll-A (RFU)`) %>%
droplevels()
rep_lin2 <- all_easylinear(Chl ~ Time|Round+Treatment+Rep, rep2, h = 7, quota = 1)
par(mfrow = c(3, 3))
par(mar = c(2, 2, 2, 1))
plot(rep_lin2, log = "y")
summary(rep_lin2)## $`2:Replete:1`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## 0.18529 -0.03780 -0.16478 -0.18114 0.17936 -0.07894 0.09801
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.37352 0.11553 37.86 2.42e-07 ***
## x 0.18726 0.01602 11.69 8.05e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1695 on 5 degrees of freedom
## Multiple R-squared: 0.9647, Adjusted R-squared: 0.9576
## F-statistic: 136.6 on 1 and 5 DF, p-value: 8.053e-05
##
##
## $`2:Replete:2`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## -0.04068 0.56895 -0.48634 -0.31279 -0.05303 0.38914 -0.06524
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.36733 0.27564 15.844 1.82e-05 ***
## x 0.20689 0.03822 5.412 0.00291 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.4045 on 5 degrees of freedom
## Multiple R-squared: 0.8542, Adjusted R-squared: 0.825
## F-statistic: 29.29 on 1 and 5 DF, p-value: 0.002913
##
##
## $`2:Replete:3`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## -0.17429 0.11485 0.12945 0.04594 -0.14406 0.10395 -0.07585
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.38836 0.22011 15.39 2.1e-05 ***
## x 0.24033 0.01335 18.01 9.7e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1412 on 5 degrees of freedom
## Multiple R-squared: 0.9848, Adjusted R-squared: 0.9818
## F-statistic: 324.3 on 1 and 5 DF, p-value: 9.701e-06
## y0 y0_lm mumax lag
## 2:Replete:1 95.47 79.32272 0.1872621 0.9894552
## 2:Replete:2 75.69 78.83269 0.2068881 -0.1966365
## 2:Replete:3 51.77 29.61731 0.2403281 2.3237062
## 2:Replete:1.r2 2:Replete:2.r2 2:Replete:3.r2
## 0.9646965 0.8542034 0.9848147
Specific Growth Rates Replete trial 3 using all_easylinear
rep3 <- fluor %>%
filter(Treatment == "Replete", Round %in% c(3)) %>%
rename(Chl = `Chlorophyll-A (RFU)`) %>%
droplevels()
rep_lin3 <- all_easylinear(Chl ~ Time|Round+Treatment+Rep, rep3, h = 7, quota = 1)
par(mfrow = c(3, 3))
par(mar = c(2, 2, 2, 1))
plot(rep_lin3, log = "y")
summary(rep_lin3)## $`3:Replete:1`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## 0.056275 0.286364 -0.388386 0.001474 -0.003777 -0.212798 0.260848
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.26798 0.31594 7.179 0.000816 ***
## x 0.26516 0.02498 10.616 0.000128 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2643 on 5 degrees of freedom
## Multiple R-squared: 0.9575, Adjusted R-squared: 0.949
## F-statistic: 112.7 on 1 and 5 DF, p-value: 0.0001282
##
##
## $`3:Replete:2`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## 0.01109 -0.07478 -0.04567 0.22757 -0.09197 -0.00871 -0.01753
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.69027 0.20383 8.292 0.000416 ***
## x 0.25682 0.01105 23.232 2.75e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.117 on 5 degrees of freedom
## Multiple R-squared: 0.9908, Adjusted R-squared: 0.989
## F-statistic: 539.7 on 1 and 5 DF, p-value: 2.75e-06
##
##
## $`3:Replete:3`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## 0.37438 -0.38671 -0.12656 -0.27670 0.47488 0.07385 -0.13314
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.09308 0.42459 4.93 0.004361 **
## x 0.30144 0.03357 8.98 0.000286 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3552 on 5 degrees of freedom
## Multiple R-squared: 0.9416, Adjusted R-squared: 0.9299
## F-statistic: 80.65 on 1 and 5 DF, p-value: 0.0002856
## y0 y0_lm mumax lag
## 3:Replete:1 46.52 9.659839 0.2651609 5.928118
## 3:Replete:2 77.27 5.420959 0.2568156 10.346075
## 3:Replete:3 61.72 8.109824 0.3014418 6.732749
## 3:Replete:1.r2 3:Replete:2.r2 3:Replete:3.r2
## 0.9575206 0.9908211 0.9416198
Specific Growth Rates Replete trial 4 using all_easylinear:
rep4 <- fluor %>%
filter(Treatment == "Replete", Round %in% c(4)) %>%
rename(Chl = `Chlorophyll-A (RFU)`) %>%
droplevels()
rep_lin4 <- all_easylinear(Chl ~ Time|Round+Treatment+Rep, rep4, h = 7, quota = 1)
par(mfrow = c(3, 3))
par(mar = c(2, 2, 2, 1))
plot(rep_lin4, log = "y")
summary(rep_lin4)## $`4:Replete:1`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## -0.09195 0.19036 0.07953 -0.20968 -0.15320 0.21726 -0.03232
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.00947 0.25121 11.98 7.15e-05 ***
## x 0.21058 0.01725 12.21 6.53e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1826 on 5 degrees of freedom
## Multiple R-squared: 0.9675, Adjusted R-squared: 0.961
## F-statistic: 149 on 1 and 5 DF, p-value: 6.529e-05
##
##
## $`4:Replete:2`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## 0.11124 0.04917 -0.34079 0.14249 -0.10230 0.22701 -0.08682
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.08172 0.28910 10.66 0.000126 ***
## x 0.20366 0.01986 10.26 0.000151 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2101 on 5 degrees of freedom
## Multiple R-squared: 0.9546, Adjusted R-squared: 0.9456
## F-statistic: 105.2 on 1 and 5 DF, p-value: 0.0001513
##
##
## $`4:Replete:3`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## 0.14336 0.08611 -0.35495 0.21552 -0.30373 0.09004 0.12366
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.07289 0.43805 4.732 0.005187 **
## x 0.17495 0.02376 7.364 0.000725 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2514 on 5 degrees of freedom
## Multiple R-squared: 0.9156, Adjusted R-squared: 0.8987
## F-statistic: 54.23 on 1 and 5 DF, p-value: 0.0007253
## y0 y0_lm mumax lag
## 4:Replete:1 78.16 20.276716 0.2105809 6.407442
## 4:Replete:2 157.37 21.795906 0.2036571 9.706892
## 4:Replete:3 87.34 7.947773 0.1749477 13.700762
## 4:Replete:1.r2 4:Replete:2.r2 4:Replete:3.r2
## 0.9675256 0.9546300 0.9155837
Combine all specific growth rate values from each treatment:
extract_results <- function(model, treatment_name) {
coefs <- coef(model)
r2vals <- rsquared(model)
df <- as.data.frame(coefs) %>%
tibble::rownames_to_column("Group") %>%
mutate(
Treatment = treatment_name,
r2 = r2vals
) %>%
# Split "Group" into Round and Rep (Group looks like "1:-Thiamine:2")
tidyr::separate(Group, into = c("Round", "Treatment_label", "Rep"), sep = ":", remove = FALSE) %>%
mutate(Round = as.integer(Round))
return(df)
}
# --- Extract Biotin rounds ---
df_bio1 <- extract_results(bio_lin1, "-Biotin")
df_bio2 <- extract_results(bio_lin2, "-Biotin")
df_bio3 <- extract_results(bio_lin3, "-Biotin")
df_bio4 <- extract_results(bio_lin4, "-Biotin")
df_bio <- bind_rows(df_bio1, df_bio2, df_bio3, df_bio4)
# --- Extract Thiamine rounds ---
df_thia1 <- extract_results(thia_lin1, "-Thiamine")
df_thia2 <- extract_results(thia_lin2, "-Thiamine")
df_thia3 <- extract_results(thia_lin3, "-Thiamine")
df_thia4 <- extract_results(thia_lin4, "-Thiamine")
df_thia <- bind_rows(df_thia1, df_thia2, df_thia3, df_thia4)
# --- Extract Replete rounds ---
df_rep1 <- extract_results(rep_lin1, "Replete")
df_rep2 <- extract_results(rep_lin2, "Replete")
df_rep3 <- extract_results(rep_lin3, "Replete")
df_rep4 <- extract_results(rep_lin4, "Replete")
df_rep <- bind_rows(df_rep1, df_rep2, df_rep3, df_rep4)
# --- Combine all treatments ---
all_lin_results <- bind_rows(df_bio, df_thia, df_rep)
# Optional: ensure Round is treated as a factor
all_lin_results$Round <- factor(all_lin_results$Round)Run an ANOVA on specific growth rate values and Tukey HSD post-hoc test for treatment:
all_lin_results$Round <- factor(all_lin_results$Round)
linear <- aov(mumax ~ Treatment , data = all_lin_results)
summary(linear)## Df Sum Sq Mean Sq F value Pr(>F)
## Treatment 2 0.00517 0.002586 1.861 0.172
## Residuals 33 0.04587 0.001390
## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = mumax ~ Treatment, data = all_lin_results)
##
## $Treatment
## diff lwr upr p adj
## -Thiamine--Biotin 0.005533672 -0.03181359 0.042880937 0.9298833
## Replete--Biotin -0.022203946 -0.05955121 0.015143319 0.3234309
## Replete--Thiamine -0.027737618 -0.06508488 0.009609647 0.1780243
Plot specific growth rates for Pyro grown in -Thiamine, -Biotin, and Replete:
all_lin_results$Treatment <- factor(
all_lin_results$Treatment,
levels = c("-Cobalamin", "-Thiamine", "-Biotin","Replete")
)
mu <- ggplot(all_lin_results, aes(x = Treatment, y = mumax, color = Treatment)) +
# Raw points
geom_point(size = 3, shape = 16, position = position_jitter(width = 0.1)) +
# Mean ± SE
stat_summary(fun = mean, geom = "point", shape = 1, size = 5, stroke = 1.2) +
stat_summary(fun.data = mean_se, geom = "errorbar", width = 0.2) +
scale_color_manual(values = c(
"-Cobalamin" = "#102E50",
"Replete" = "#F5C45E",
"-Thiamine" = "deeppink4",
"-Biotin" = "#667028"
)) +
scale_x_discrete(labels = c(
"-Cobalamin" = expression("\u2212B"[12]),
"Replete" = "Replete",
"-Thiamine" = expression("\u2212Thiamine"),
"-Biotin" = expression("\u2212Biotin")
)) +
theme_minimal() +
theme(
legend.position = "none",
panel.grid = element_blank(),
panel.border = element_rect(fill = NA, color = "black"),
axis.title = element_text( size = 20),
axis.text = element_text( size = 20, color = "black"),
axis.text.x = element_text(angle = 45, hjust = 1),
strip.text = element_text( size = 20, color = "black"),
legend.title = element_blank(),
legend.text = element_blank(),
legend.key = element_blank()
) +
labs(
x = "Media Treatment",
y = expression(paste("(", mu, ") RFU" * "\u00B7" * "day"^-1, " "))
) + guides(color = "none", fill = "none", shape = "none")
mubo <-
LD_OTB_all_growth /
(max_bio_LD +
mu ) +
plot_layout(guides = "collect") &
theme(
legend.position = "bottom",
panel.spacing = unit(.1, "lines"),
plot.margin = margin(15, 15, 15, 15),
axis.title.y = element_text(margin = margin(r = 12)),
axis.title.x = element_text(margin = margin(t = 8)),
legend.margin = margin(t = 4, b = 4),
legend.box.margin = margin(5, 5, 5, 5),
axis.text.x = element_text(size = 16),
)
boChlorophyll-a fluorescence (RFU) data from P. bahamense cultures were used to visualize growth patterns, quantify maximum biomass, and estimate specific growth rates under media conditions deficient in media lacking one targeted B vitamin. Specific growth rates were computed using the “alleasy_linear” function of the “growthrates” package in R Studio (Hall et al., 2014; RStudio Team, 2020).
Import and Parse the Dataset:
b12 <-as.data.frame(read_csv("b12_levels.csv"))
# Remove periods and convert to numeric
b12$Time <- as.numeric(sub("T", "", b12$Time))
# Convert Date column to Date type
b12$Date <- mdy(b12$Date)
b12$Treatment <- unlist(b12$Treatment)
#covert timepoints to hours
b12$Time <- b12$Time * 48 # hours
b12$Time <- b12$Time / 24
# put treatment level in desired order
b12 <- b12 %>%
mutate(Treatment = factor(Treatment, levels = c(
"369pM", # Place the levels in desired order
"36pM",
"3.69pM",
"0.369pM",
"0.0369pM"
)))Prepare data for plotting:
levels_order <- c("0.0369pM", "0.369pM", "3.69pM", "36pM", "369pM")
b12$Treatment <- factor(b12$Treatment, levels = levels_order)
# Named color palette
pp_palette <- c(
"0.0369pM" = "lightblue",
"0.369pM" = "#99c2e6",
"3.69pM" = "#66a3d2",
"36pM" = "#3380bf",
"369pM" = "#004080"
)growth <- ggplot(b12, aes(x = Time, y = `Chlorophyll-A (RFU)`, color = Treatment)) +
geom_point() +
geom_smooth(aes(color = Treatment), se = FALSE) +
facet_grid(Round ~ Treatment) +
labs(
x = "Time (days)",
y = expression(paste("Chlorophyll-", italic("a"), " fluorescence (RFU)"))
) +
scale_y_continuous(labels = scales::comma) +
scale_color_manual(values = pp_palette) +
theme_minimal() +
theme(
panel.border = element_rect(fill = NA, color = "black"),
panel.grid = element_blank(),
axis.text.y = element_text(size = 20, color = "black"),
axis.text.x = element_text(size = 12, color = "black"),
axis.title.y = element_text(size = 20, margin = margin(t = 10)),
axis.title.x = element_text(size = 20, margin = margin(t = 10)),
strip.text.x = element_text(size = 16),
strip.text = element_text( size = 16, color = "black"),
legend.title = element_text( size = 16, color = "black"),
legend.text = element_text( size = 16, color = "black"),
legend.position = "bottom",
axis.ticks = element_line(color = "black", linewidth = 0.5),
axis.ticks.length = unit(4, "pt"),
panel.spacing = unit(.15, "cm")
)
growthmax_bio <- b12 %>%
filter(Treatment %in% c(Treatment, levels = c("369pM", "36pM", "3.69pM", "0.369pM", "0.0369pM"))) %>%
group_by(Treatment, Bio_Rep) %>%
summarise(Max_Chl = max(`Chlorophyll-A (RFU)`, na.rm = TRUE)) %>%
ungroup()
print(max_bio)## # A tibble: 15 × 3
## Treatment Bio_Rep Max_Chl
## <fct> <dbl> <dbl>
## 1 0.0369pM 1 842.
## 2 0.0369pM 2 541.
## 3 0.0369pM 3 928.
## 4 0.369pM 1 1968.
## 5 0.369pM 2 1897.
## 6 0.369pM 3 1866.
## 7 3.69pM 1 1923.
## 8 3.69pM 2 3929.
## 9 3.69pM 3 2146.
## 10 36pM 1 5580.
## 11 36pM 2 6041.
## 12 36pM 3 4801.
## 13 369pM 1 5287.
## 14 369pM 2 4424.
## 15 369pM 3 5539.
Run an ANOVA followed by a Tukey HSD test for the treatment factor
## Df Sum Sq Mean Sq F value Pr(>F)
## Treatment 4 49705730 12426432 31.27 1.25e-05 ***
## Residuals 10 3973244 397324
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = Max_Chl ~ Treatment, data = max_bio)
##
## $Treatment
## diff lwr upr p adj
## 0.369pM-0.0369pM 1139.8033 -554.0107 2833.617 0.2494204
## 3.69pM-0.0369pM 1895.5967 201.7826 3589.411 0.0272120
## 36pM-0.0369pM 4703.7467 3009.9326 6397.561 0.0000277
## 369pM-0.0369pM 4312.9433 2619.1293 6006.757 0.0000597
## 3.69pM-0.369pM 755.7933 -938.0207 2449.607 0.6023509
## 36pM-0.369pM 3563.9433 1870.1293 5257.757 0.0003042
## 369pM-0.369pM 3173.1400 1479.3259 4866.954 0.0007806
## 36pM-3.69pM 2808.1500 1114.3359 4501.964 0.0020036
## 369pM-3.69pM 2417.3467 723.5326 4111.161 0.0058828
## 369pM-36pM -390.8033 -2084.6174 1303.011 0.9366508
Plot maximum biomass:
# Define factor levels (low to high B12)
levels_order <- c("0.0369pM", "0.369pM", "3.69pM", "36pM", "369pM")
max_bio$Treatment <- factor(trimws(max_bio$Treatment), levels = levels_order)
# Calculate summary statistics
summary_stats <- max_bio %>%
group_by(Treatment) %>%
summarise(
mean_chl = mean(Max_Chl, na.rm = TRUE),
sd_chl = sd(Max_Chl, na.rm = TRUE),
.groups = "drop"
)
summary_stats$Treatment <- factor(summary_stats$Treatment, levels = levels_order)
# Define color palette
pp_palette <- c(
"0.0369pM" = "lightblue",
"0.369pM" = "#99c2e6",
"3.69pM" = "#66a3d2",
"36pM" = "#3380bf",
"369pM" = "#004080"
)
max_bio_plot <- ggplot(summary_stats, aes(x = Treatment)) +
# Replicate points
geom_point(
data = max_bio,
aes(y = Max_Chl, color = Treatment),
size = 3,
alpha = 1,
position = position_jitter(width = 0.1)
) +
stat_summary(
data = max_bio,
aes(y = Max_Chl, color = Treatment),
fun = mean,
geom = "point",
shape = 1,
size = 5,
stroke = 1.2
) +
stat_summary(
data = max_bio,
aes(y = Max_Chl, color = Treatment),
fun.data = mean_sdl,
fun.args = list(mult = 1),
geom = "errorbar",
width = 0.2
)+
# Keep your palette but remove the scale (so no legend)
scale_color_manual(values = pp_palette, guide = "none") +
scale_y_continuous(labels = scales::comma) +
coord_cartesian(ylim = c(0, 10000)) +
ylab(expression(paste("Maximum RFU"))) +
xlab(expression(B[12]~Concentration~"(pM)")) +
theme_minimal() +
theme(
panel.border = element_rect(fill = NA, color = "black"),
panel.grid = element_blank(),
axis.text.y = element_text(size = 20, color = "black"),
axis.text.x = element_text(size = 20, color = "black", angle = 45, hjust = 1),
axis.title.y = element_text(size = 20),
axis.title.x = element_text(size = 20),
legend.position = "none",
axis.ticks = element_line(color = "black", linewidth = 0.5),
axis.ticks.length = unit(4, "pt")
)
max_bio_plotRemove outliers:
# ️⃣ Compute T0 baseline per replicate, treatment, and round
baseline <- aggregate(
`Chlorophyll-A (RFU)` ~ Bio_Rep + Treatment + Round,
data = subset(b12, Time == 0),
FUN = mean,
na.rm = TRUE
)
names(baseline)[names(baseline) == "Chlorophyll-A (RFU)"] <- "Chl_T0"
# Merge baseline back into the full dataset
b12_merged <- merge(b12, baseline, by = c("Bio_Rep", "Treatment", "Round"), all.x = TRUE)
# Keep all T0 points, and remove any post-T0 points where replicate < its T0
b12_clean <- subset(
b12_merged,
Time == 0 | `Chlorophyll-A (RFU)` >= Chl_T0
)Specific growth rate with all easy linear for 0.0369PM:
pM_1 <- b12_clean %>%
filter(
Treatment %in% c("0.0369pM"),
Round %in% c(1)
) %>%
rename(Chl = `Chlorophyll-A (RFU)`) %>%
arrange(Time, Bio_Rep)
mu1 <- all_easylinear(Chl ~ Time | Round + Treatment + Bio_Rep, pM_1, h =4, quota = 1)
par(mfrow = c(3, 3))
par(mar = c(2, 2, 2, 1))
plot(mu1, log = "y")
summary(mu1)## $`1:0.0369pM:1`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4
## 0.21354 -0.37851 0.11642 0.04856
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.72824 0.52579 8.993 0.0121 *
## x 0.19590 0.07155 2.738 0.1115
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.32 on 2 degrees of freedom
## Multiple R-squared: 0.7894, Adjusted R-squared: 0.6841
## F-statistic: 7.496 on 1 and 2 DF, p-value: 0.1115
##
##
## $`1:0.0369pM:2`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4
## -0.10998 0.29391 -0.05595 -0.12799
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 5.39799 0.22595 23.890 0.00175 **
## x 0.07310 0.02382 3.069 0.09177 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2429 on 2 degrees of freedom
## Multiple R-squared: 0.8249, Adjusted R-squared: 0.7373
## F-statistic: 9.421 on 1 and 2 DF, p-value: 0.09177
##
##
## $`1:0.0369pM:3`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4
## 0.29094 -0.55228 0.23175 0.02959
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.5299 0.5772 7.848 0.0159 *
## x 0.2125 0.1054 2.017 0.1813
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.4713 on 2 degrees of freedom
## Multiple R-squared: 0.6703, Adjusted R-squared: 0.5055
## F-statistic: 4.067 on 1 and 2 DF, p-value: 0.1813
## y0 y0_lm mumax lag
## 1:0.0369pM:1 145.02 113.09593 0.19590188 1.2691828
## 1:0.0369pM:2 197.95 220.96199 0.07310236 -1.5044141
## 1:0.0369pM:3 97.44 92.75275 0.21252290 0.2319725
## 1:0.0369pM:1.r2 1:0.0369pM:2.r2 1:0.0369pM:3.r2
## 0.7893916 0.8248814 0.6703408
Specific growth rate with all easy linear for 0.369pM:
pM_2 <- b12_clean %>%
filter(Treatment %in% c("0.369pM"),
Round %in% c(1),
Time >= 0) %>% # <-- exclude T0
rename(Chl = `Chlorophyll-A (RFU)`) %>%
arrange(Time, Bio_Rep)
mu2 <-all_easylinear(Chl ~ Time|Treatment+Bio_Rep, pM_2, h = 4, quota = 1)
par(mfrow = c(3, 3))
par(mar = c(2, 2, 2, 1))
plot(mu2, log = "y") # line only
summary(mu2)## $`0.369pM:1`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4
## -0.03871 0.14944 -0.14406 0.03332
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.72896 0.22995 20.566 0.00236 **
## x 0.23253 0.02555 9.101 0.01186 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1512 on 2 degrees of freedom
## Multiple R-squared: 0.9764, Adjusted R-squared: 0.9646
## F-statistic: 82.83 on 1 and 2 DF, p-value: 0.01186
##
##
## $`0.369pM:2`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4
## -0.03225 -0.15070 0.39814 -0.21519
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 5.12741 0.70102 7.314 0.0182 *
## x 0.20224 0.07559 2.675 0.1159
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3381 on 2 degrees of freedom
## Multiple R-squared: 0.7816, Adjusted R-squared: 0.6724
## F-statistic: 7.157 on 1 and 2 DF, p-value: 0.1159
##
##
## $`0.369pM:3`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4
## -0.26163 0.40004 -0.01518 -0.12322
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.21953 0.72414 5.827 0.0282 *
## x 0.28628 0.07809 3.666 0.0670 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3492 on 2 degrees of freedom
## Multiple R-squared: 0.8705, Adjusted R-squared: 0.8057
## F-statistic: 13.44 on 1 and 2 DF, p-value: 0.06701
## y0 y0_lm mumax lag
## 0.369pM:1 248.73 113.17754 0.2325283 3.386298
## 0.369pM:2 375.27 168.57963 0.2022351 3.956966
## 0.369pM:3 182.06 68.00143 0.2862760 3.440063
## 0.369pM:1.r2 0.369pM:2.r2 0.369pM:3.r2
## 0.9764234 0.7815961 0.8704739
Specific growth rate with all easy linear for 3.69PM:
pM_3 <- b12_clean %>%
filter(Treatment %in% c("3.69pM"),
Round %in% c(1),
Time >= 0) %>%
rename(Chl = `Chlorophyll-A (RFU)`) %>%
arrange(Time, Bio_Rep)
mu3 <-all_easylinear(Chl ~ Time|Round+Treatment+Bio_Rep, pM_3, h =5, quota = 1)
par(mfrow = c(3, 3))
par(mar = c(2, 2, 2, 1))
plot(mu3 , log = "y")
summary(mu3)## $`1:3.69pM:1`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5
## -0.03126 -0.21999 0.32699 0.13101 -0.20676
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.60153 0.36017 12.776 0.00103 **
## x 0.24850 0.04245 5.854 0.00994 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2685 on 3 degrees of freedom
## Multiple R-squared: 0.9195, Adjusted R-squared: 0.8927
## F-statistic: 34.28 on 1 and 3 DF, p-value: 0.009935
##
##
## $`1:3.69pM:2`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5
## 0.100883 -0.081917 -0.005946 -0.145891 0.132871
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 5.40262 0.10570 51.11 1.65e-05 ***
## x 0.21586 0.02158 10.01 0.00213 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1365 on 3 degrees of freedom
## Multiple R-squared: 0.9709, Adjusted R-squared: 0.9612
## F-statistic: 100.1 on 1 and 3 DF, p-value: 0.002126
##
##
## $`1:3.69pM:3`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5
## -0.076428 0.082299 0.033975 -0.009133 -0.030712
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.8058 0.1154 41.65 3.05e-05 ***
## x 0.2069 0.0111 18.63 0.000337 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.07023 on 3 degrees of freedom
## Multiple R-squared: 0.9914, Adjusted R-squared: 0.9886
## F-statistic: 347.1 on 1 and 3 DF, p-value: 0.0003375
## y0 y0_lm mumax lag
## 1:3.69pM:1 142.16 99.63649 0.2485048 1.4302532
## 1:3.69pM:2 245.55 221.98664 0.2158566 0.4673628
## 1:3.69pM:3 271.24 122.21605 0.2068770 3.8535635
## 1:3.69pM:1.r2 1:3.69pM:2.r2 1:3.69pM:3.r2
## 0.9195173 0.9708995 0.9914315
Specific growth rate with all easy linear for 36PM:
pM_4 <- b12_clean %>%
filter(Treatment %in% c("36pM"),
Round %in% c(1),
Time >= 0) %>% # <-- exclude T0
rename(Chl = `Chlorophyll-A (RFU)`) %>%
arrange(Time, Bio_Rep)
mu4 <-all_easylinear(Chl ~ Time|Round+Treatment+Bio_Rep, pM_4, h = 5, quota = 1)
par(mfrow = c(3, 3))
par(mar = c(2, 2, 2, 1))
plot(mu4 , log = "y")
summary(mu4)## $`1:36pM:1`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5
## 0.02088 -0.04810 0.12932 -0.19786 0.09576
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.69402 0.24700 19.00 0.000318 ***
## x 0.23918 0.02377 10.06 0.002090 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1503 on 3 degrees of freedom
## Multiple R-squared: 0.9712, Adjusted R-squared: 0.9616
## F-statistic: 101.3 on 1 and 3 DF, p-value: 0.00209
##
##
## $`1:36pM:2`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5
## -0.5827 0.4830 0.2809 0.1390 -0.3201
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.92981 0.81583 6.043 0.00909 **
## x 0.18334 0.06122 2.995 0.05791 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.5078 on 3 degrees of freedom
## Multiple R-squared: 0.7494, Adjusted R-squared: 0.6658
## F-statistic: 8.969 on 1 and 3 DF, p-value: 0.05791
##
##
## $`1:36pM:3`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5
## 0.05819 0.08475 -0.32792 0.16882 0.01616
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.22326 0.36356 11.617 0.00137 **
## x 0.27122 0.03498 7.753 0.00446 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2213 on 3 degrees of freedom
## Multiple R-squared: 0.9525, Adjusted R-squared: 0.9366
## F-statistic: 60.11 on 1 and 3 DF, p-value: 0.004463
## y0 y0_lm mumax lag
## 1:36pM:1 330.46 109.29127 0.2391784 4.626125
## 1:36pM:2 205.58 138.35322 0.1833369 2.160096
## 1:36pM:3 104.26 68.25567 0.2712184 1.561941
## 1:36pM:1.r2 1:36pM:2.r2 1:36pM:3.r2
## 0.9712271 0.7493526 0.9524610
Specific growth rate with all easy linear for 369PM:
pM_5 <- b12_clean %>%
filter(Treatment %in% c("369pM"),
Round %in% c(1),
Time >= 0) %>%
rename(Chl = `Chlorophyll-A (RFU)`) %>%
arrange(Time, Bio_Rep)
mu5 <-all_easylinear(Chl ~ Time|Round+Treatment+Bio_Rep, pM_5, h = 5, quota = 1)
par(mfrow = c(3, 3))
par(mar = c(2, 2, 2, 1))
plot(mu5 , log = "y")
summary(mu5)## $`1:369pM:1`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5
## -0.08633 0.17016 -0.04861 0.01842 -0.05363
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.74188 0.15680 30.24 7.94e-05 ***
## x 0.25385 0.01538 16.51 0.000484 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1183 on 3 degrees of freedom
## Multiple R-squared: 0.9891, Adjusted R-squared: 0.9855
## F-statistic: 272.6 on 1 and 3 DF, p-value: 0.0004836
##
##
## $`1:369pM:2`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5
## 0.221300 -0.358207 0.001927 -0.035733 0.170713
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.6394 0.3488 13.303 0.000918 ***
## x 0.2400 0.0342 7.017 0.005946 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2631 on 3 degrees of freedom
## Multiple R-squared: 0.9426, Adjusted R-squared: 0.9234
## F-statistic: 49.23 on 1 and 3 DF, p-value: 0.005946
##
##
## $`1:369pM:3`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5
## -0.26764 0.51741 -0.14472 -0.19223 0.08718
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.27559 0.49285 8.675 0.00322 **
## x 0.25738 0.05808 4.431 0.02135 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3673 on 3 degrees of freedom
## Multiple R-squared: 0.8675, Adjusted R-squared: 0.8233
## F-statistic: 19.64 on 1 and 3 DF, p-value: 0.02135
## y0 y0_lm mumax lag
## 1:369pM:1 286.00 114.64947 0.2538508 3.600983
## 1:369pM:2 321.18 103.48174 0.2399582 4.720016
## 1:369pM:3 106.96 71.92225 0.2573820 1.541947
## 1:369pM:1.r2 1:369pM:2.r2 1:369pM:3.r2
## 0.9891142 0.9425652 0.8674686
Make dataframe with fits:
make_df <- function(model, round_num, experiment_name) {
coefs <- coef(model)
if (is.matrix(coefs)) {
df <- as.data.frame(coefs) %>%
tibble::rownames_to_column("Group") %>%
mutate(
Round = round_num,
Experiment = experiment_name,
r2 = rsquared(model)
)
} else {
df <- data.frame(
Group = names(coefs),
Estimate = as.numeric(coefs),
Round = round_num,
Experiment = experiment_name,
r2 = rep(rsquared(model), length(coefs))
)
}
return(df)
}
# Example usage for your single fits
df1 <- make_df(mu1, 1, "0.0369pM")
df2 <- make_df(mu2, 1, "0.369pM")
df3 <- make_df(mu3, 1, "3.69pM")
df4 <- make_df(mu4, 1, "36pM")
df5 <- make_df(mu5, 1, "369pM")
# Combine all results
all_results <- bind_rows(df1, df2, df3, df4, df5)
# Check
all_results## Group y0 y0_lm mumax lag Round Experiment
## 1 1:0.0369pM:1 145.02 113.09593 0.19590188 1.2691828 1 0.0369pM
## 2 1:0.0369pM:2 197.95 220.96199 0.07310236 -1.5044141 1 0.0369pM
## 3 1:0.0369pM:3 97.44 92.75275 0.21252290 0.2319725 1 0.0369pM
## 4 0.369pM:1 248.73 113.17754 0.23252831 3.3862982 1 0.369pM
## 5 0.369pM:2 375.27 168.57963 0.20223514 3.9569659 1 0.369pM
## 6 0.369pM:3 182.06 68.00143 0.28627601 3.4400634 1 0.369pM
## 7 1:3.69pM:1 142.16 99.63649 0.24850476 1.4302532 1 3.69pM
## 8 1:3.69pM:2 245.55 221.98664 0.21585664 0.4673628 1 3.69pM
## 9 1:3.69pM:3 271.24 122.21605 0.20687699 3.8535635 1 3.69pM
## 10 1:36pM:1 330.46 109.29127 0.23917836 4.6261254 1 36pM
## 11 1:36pM:2 205.58 138.35322 0.18333686 2.1600963 1 36pM
## 12 1:36pM:3 104.26 68.25567 0.27121840 1.5619415 1 36pM
## 13 1:369pM:1 286.00 114.64947 0.25385084 3.6009826 1 369pM
## 14 1:369pM:2 321.18 103.48174 0.23995820 4.7200161 1 369pM
## 15 1:369pM:3 106.96 71.92225 0.25738195 1.5419469 1 369pM
## r2
## 1 0.7893916
## 2 0.8248814
## 3 0.6703408
## 4 0.9764234
## 5 0.7815961
## 6 0.8704739
## 7 0.9195173
## 8 0.9708995
## 9 0.9914315
## 10 0.9712271
## 11 0.7493526
## 12 0.9524610
## 13 0.9891142
## 14 0.9425652
## 15 0.8674686
ANOVA and post-hoc multiple comparisons of mumax among varying media treatments:
## Df Sum Sq Mean Sq F value Pr(>F)
## Experiment 4 0.01503 0.003757 1.85 0.196
## Residuals 10 0.02031 0.002031
## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = mumax ~ Experiment, data = all_results)
##
## $Experiment
## diff lwr upr p adj
## 0.369pM-0.0369pM 0.079837439 -0.04125361 0.2009285 0.2651098
## 3.69pM-0.0369pM 0.063237082 -0.05785396 0.1843281 0.4652490
## 369pM-0.0369pM 0.089887949 -0.03120310 0.2109790 0.1807721
## 36pM-0.0369pM 0.070735492 -0.05035555 0.1918265 0.3656496
## 3.69pM-0.369pM -0.016600357 -0.13769140 0.1044907 0.9901101
## 369pM-0.369pM 0.010050511 -0.11104054 0.1311416 0.9985578
## 36pM-0.369pM -0.009101947 -0.13019299 0.1119891 0.9990215
## 369pM-3.69pM 0.026650867 -0.09444018 0.1477419 0.9459297
## 36pM-3.69pM 0.007498410 -0.11359264 0.1285895 0.9995435
## 36pM-369pM -0.019152457 -0.14024350 0.1019386 0.9831961
T-test between 0.0369pM and all other concentrations:
all_results$Group <- ifelse(all_results$Experiment == "0.0369pM",
"A",
"B")
all_results$Group <- factor(all_results$Group)
t.test(mumax~ Group, data = all_results)##
## Welch Two Sample t-test
##
## data: mumax by Group
## t = -1.6943, df = 2.1577, p-value = 0.2231
## alternative hypothesis: true difference in means between group A and group B is not equal to 0
## 95 percent confidence interval:
## -0.2558391 0.1039901
## sample estimates:
## mean in group A mean in group B
## 0.1605090 0.2364335
Plot mumax with all easy linear data:
# Define the order of treatments from least to greatest
levels_order <- c("0.0369pM", "0.369pM", "3.69pM", "36pM", "369pM")
# Convert Treatment column to a factor with this order
all_results$Experiment <- factor(all_results$Experiment, levels = levels_order)
# calculate standard deviation and mean
b12_summary <- all_results %>%
group_by(Experiment) %>%
summarise(
mean_mumax = mean(mumax, na.rm = TRUE),
se_mumax = sd(mumax, na.rm = TRUE)/sqrt(n()),
.groups = "drop"
)
mu_b12 <- ggplot(all_results, aes(x = Experiment, y = mumax, color = Experiment)) +
# Raw replicate points (jittered)
geom_point(size = 3, shape = 16,
position = position_jitter(width = 0.1), alpha = 1) +
# Mean per Experiment+Round (open circle)
stat_summary(fun = mean, geom = "point", shape = 1, size = 5, stroke = 1.2) +
# Error bars (standard error)
stat_summary(fun.data = mean_se, geom = "errorbar",
width = 0.2) +
# Color palette — WITH legend removed
scale_color_manual(
values = c(
"0.0369pM" = "lightblue",
"0.369pM" = "#99c2e6",
"3.69pM" = "#66a3d2",
"36pM" = "#3380bf",
"369pM" = "#004080"
),
guide = "none" # ← THIS removes the legend safely
) +
# Safe visible range
coord_cartesian(ylim = c(0, 0.4)) +
# Labels + theme
labs(
x = expression(B[12]~Concentration~"(pM)"),
y = expression(paste(mu, " (RFU" %*% "day"^-1, ")"))
) +
theme_minimal() +
theme(
legend.position = "none",
panel.grid = element_blank(),
panel.border = element_rect(fill = NA, color = "black"),
axis.title = element_text(size = 20),
axis.text = element_text(size = 20, color = "black"),
axis.text.x = element_text(size = 20, color = "black", angle = 45, hjust = 1),
strip.text = element_text(size = 20, color = "black"),
axis.ticks = element_line(color = "black", linewidth = 0.5),
axis.ticks.length = unit(4, "pt")
)
mu_b12# Ensure BOTH subplots have no legend inside patchwork
max_bio_nl <- max_bio_plot + theme(legend.position = "none")
mu_nl <- mu_b12 + theme(legend.position = "none")
# Build the layout
bo <-
(
growth /
(max_bio_nl | mu_nl)
) +
plot_layout(
widths = c(4, 2), # your original width ratio
guides = "collect" # ← THIS is what actually gathers one legend
) &
theme(
legend.position = "bottom",
legend.box = "horizontal",
panel.spacing = unit(0.4, "lines"),
plot.margin = margin(15, 15, 15, 15),
axis.title.y = element_text(margin = margin(r = 12)),
axis.title.x = element_text(margin = margin(t = 8)),
legend.margin = margin(t = 4, b = 4),
legend.box.margin = margin(5, 5, 5, 5),
legend.ticks = element_blank(),
axis.text.x = element_text(size = 16)
)
boLow-diveristy P. bahamense culture was co-cultured with bacteria recruited from OTB (Old Tampa Bay) water, Indian River Lagoon (IRL) culture, and IRL surface seawater and grown in replete media and media missing B12
Import data:
OTB water:
OTB <- as.data.frame(read_csv("OTB_1006.csv"))
OTB$Time <- as.numeric(sub('.', '', OTB$Time))
OTB$Date <- ymd(OTB$Date)
OTB$Time <- OTB$Time * 48
OTB %<>%
filter(Treatment %in% c("-Cobalamin", "Replete")) %>%
mutate(Treatment = factor(Treatment, levels = c("-Cobalamin", "Replete")))
OTB$Time_in_Days <- OTB$Time / 24IRL water:
IRLrec<-as.data.frame(read_csv("IRL_1006.csv"))
IRLrec$Time <- as.numeric(sub('.', '', IRLrec$Time))
IRLrec$Date <- ymd(IRLrec$Date)
IRLrec$Time <- IRLrec$Time * 48
IRLrec$Time_in_Days <- IRLrec$Time / 24 # Convert Time from hours to daysIRL culture:
IRLcul <- as.data.frame(read_csv("1003_1006.csv"))
IRLcul$Time <- as.numeric(gsub("\\.", "", IRLcul$Time)) # Remove all periods
IRLcul$Date <- ymd(IRLcul$Date)
IRLcul$Time <- IRLcul$Time * 48
IRLcul$Time_in_Days <- IRLcul$Time / 24 Combine all datasets into one with a new ‘Dataset’ column:
all_experiments <- bind_rows(
fluor %>% mutate(Dataset = "Low Diversity OTB"),
IRLcul %>% mutate(Dataset = "IRL culture microbiome"),
OTB %>% mutate(Dataset = "OTB water microbiome"),
IRLrec %>% mutate(Dataset = "IRL water microbiome")
) %>%mutate(
Treatment = recode(Treatment,
"-Cobalamin" = "\u2013B12",
"-Biotin" = "\u2013Biotin",
"-Thiamine" = "\u2013Thiamine"),
Treatment = factor(Treatment, levels = c("\u2013B12", "Replete"))
) %>%
filter(Round %in% c(1, 2), Treatment %in% c("\u2013B12", "Replete"))# Filter out "Low Diversity OTB"
all_experiments_filtered <- all_experiments %>%
filter(Dataset != "Low Diversity OTB") %>%
mutate(
Treatment = recode(Treatment,
"–B12" = "—B12", # convert en dash to em dash
"Replete" = "Replete"),
Treatment = factor(Treatment, levels = c("—B12", "Replete")),
Dataset = factor(Dataset, levels = c(
"IRL culture microbiome",
"OTB water microbiome",
"IRL water microbiome"
))
)
growth <- ggplot(all_experiments_filtered, aes(x = Time_in_Days, y = `Chlorophyll-A (RFU)`, color = Treatment)) +
geom_point(size = 3) +
geom_smooth(se = FALSE) +
scale_color_manual(values = c("—B12" = "#102E50", "Replete" = "#F5C45E")) +
facet_grid(Round ~ Dataset) +
labs(
x = "Time (days)",
y = expression(paste("Chlorophyll-", italic("a"), " fluorescence (RFU)"))
) +
scale_y_continuous(labels = label_comma()) +
theme_minimal() +
theme(
panel.border = element_rect(fill = NA, color = "black"),
panel.grid = element_blank(),
axis.text.y = element_text(size = 16, color = "black"),
axis.text.x = element_text(size = 12, color = "black"),
axis.title.y = element_text(size = 20, margin = margin(t = 10)),
axis.title.x = element_text(size = 20, margin = margin(t = 10)),
strip.text.x = element_text(size = 20),
strip.text = element_text(size = 20, color = "black"),
legend.title = element_text(size = 16, color = "black"),
legend.text = element_text(size = 16, color = "black"),
legend.position = "bottom",
axis.ticks = element_line(color = "black", linewidth = 0.5),
axis.ticks.length = unit(4, "pt"),
panel.spacing = unit(.1, "cm") ) +
coord_cartesian(ylim = c(0, 5000))
plot(growth)Create new variable for experimemt and combine experimental data into dataframe:
max_bio <- all_exp %>%
group_by(Treatment, Bio_Rep,Experiment, Round) %>%
summarise(Max_Chl = max(`Chlorophyll-A (RFU)`, na.rm = TRUE), .groups = "drop") %>%
arrange(Treatment, Bio_Rep)
print(max_bio)## # A tibble: 44 × 5
## Treatment Bio_Rep Experiment Round Max_Chl
## <chr> <dbl> <chr> <dbl> <dbl>
## 1 -Cobalamin 1 IRLcul 1 103.
## 2 -Cobalamin 1 IRLcul 2 136.
## 3 -Cobalamin 1 IRLwat 1 893.
## 4 -Cobalamin 1 IRLwat 2 904.
## 5 -Cobalamin 1 OTB 1 341.
## 6 -Cobalamin 1 OTB 2 87.8
## 7 -Cobalamin 2 IRLcul 1 147.
## 8 -Cobalamin 2 IRLcul 2 173.
## 9 -Cobalamin 2 IRLwat 1 674.
## 10 -Cobalamin 2 IRLwat 2 1126.
## # ℹ 34 more rows
Plot maximum biomass
# Define experiment order
experiment_order <- c("fluor", "IRLcul", "OTB", "IRLwat")
# Filter and prepare data for -Cobalamin only
plot_data <- max_bio %>%
filter(Treatment == "-Cobalamin") %>%
mutate(Experiment = factor(Experiment, levels = experiment_order))
#
# plot_summary <- summary_stats %>%
# filter(Treatment == "-Cobalamin") %>%
# mutate(Experiment = factor(Experiment, levels = experiment_order))
# Build plot
max_bio_plot <- ggplot(plot_data, aes(x = Experiment, y = Max_Chl, color = Treatment)) +
# Raw replicate points (solid)
geom_point(size = 5, shape = 16, alpha = 1,
position = position_jitter(width = 0.15)) +
stat_summary(
fun = mean,
geom = "point",
shape = 1,
size = 5,
stroke = 1.2
) +
# Error bars: mean ± SD (NOT SE)
stat_summary(
fun.data = mean_sdl,
fun.args = list(mult = 1), # 1 SD above/below mean
geom = "errorbar",
width = 0.2
) +
# Legend with open circle
scale_color_manual(
values = c("-Cobalamin" = "#102E50"),
labels = c("-Cobalamin" = expression("—B"[12])),
name = "Treatment",
guide = guide_legend(
override.aes = list(
shape = 21,
fill = "white", # open circle in legend
stroke = 1.2,
size = 5
)
)
) +
# Axis labels
labs(
x = "Source",
y = expression(paste("Maximum RFU"))
) +
# X-axis labels
scale_x_discrete(labels = c(
fluor = "Low Diversity",
IRLcul = "IRL Culture",
OTB = "OTB Seawater",
IRLwat = "IRL Seawater"
)) +
# Y-axis formatting
scale_y_continuous(labels = label_comma()) +
# # Facet by Round
# facet_grid(~Round) +
# Theme
theme_minimal() +
theme(
panel.border = element_rect(fill = NA, color = "black"),
panel.grid = element_blank(),
axis.text.y = element_text(size = 20, color = "black"),
axis.text.x = element_text(angle = 45, hjust = 1, size = 14, color = "black"),
axis.title.y = element_text(size = 20),
axis.title.x = element_text(size = 20, margin = margin(t = 10)),
strip.text.x = element_text(size = 20),
strip.text = element_text(size = 20, color = "black"),
legend.title = element_text(size = 16, color = "black"),
legend.text = element_text(size = 16, color = "black"),
legend.position = "bottom",
legend.key.size = unit(1.5, "cm"),
axis.ticks = element_line(color = "black", linewidth = 0.5),
axis.ticks.length = unit(4, "pt")
) +
# Set visible y-axis range safely
coord_cartesian(ylim = c(50, 1500))
# Print plot
max_bio_plotRun an ANOVA followed by a Tukey HSD test for the experiment factor:
max_B12 <- all_exp %>%
group_by(Experiment ,Treatment, Round, Bio_Rep) %>%
summarise(Max_Chl = max(`Chlorophyll-A (RFU)`, na.rm = TRUE)) %>%
ungroup() %>% filter(Treatment %in% c("-Cobalamin"))
aov_biomass <- aov(Max_Chl ~ Experiment, data =max_B12)
summary(aov_biomass)## Df Sum Sq Mean Sq F value Pr(>F)
## Experiment 3 2481763 827254 38.95 4.45e-08 ***
## Residuals 18 382327 21240
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = Max_Chl ~ Experiment, data = max_B12)
##
## $Experiment
## diff lwr upr p adj
## IRLcul-fluor -148.55167 -414.4357 117.3324 0.4145013
## IRLwat-fluor 679.19000 413.3060 945.0740 0.0000057
## OTB-fluor -23.03833 -288.9224 242.8457 0.9946502
## IRLwat-IRLcul 827.74167 589.9278 1065.5556 0.0000001
## OTB-IRLcul 125.51333 -112.3006 363.3272 0.4626003
## OTB-IRLwat -702.22833 -940.0422 -464.4144 0.0000007
Remove outliers from IRL culture data: outliers are considered values the drop below T0 RFU values
# Compute T0 baseline per replicate, treatment, and round
baseline <- aggregate(
`Chlorophyll-A (RFU)` ~ Bio_Rep + Treatment + Round,
data = subset(IRLcul, Time == 0),
FUN = mean,
na.rm = TRUE
)
names(baseline)[names(baseline) == "Chlorophyll-A (RFU)"] <- "Chl_T0"
# Merge baseline back into the full dataset
IRLcul_merged <- merge(IRLcul, baseline, by = c("Bio_Rep", "Treatment", "Round"), all.x = TRUE)
# Keep all T0 points, and remove any post-T0 points where replicate < its T0
IRLcul_clean <- subset(
IRLcul_merged,
Time == 0 | `Chlorophyll-A (RFU)` >= Chl_T0
)IRL culture trial 1 Replete:
IRLcul_mumax <- IRLcul_clean %>%
filter(Treatment %in% c("Replete"),
Round %in% c(1),
Time >= 0) %>% # <-- exclude T0
rename(Chl = `Chlorophyll-A (RFU)`) %>%
arrange(Time, Bio_Rep)
IRLcul1 <- all_easylinear(Chl ~ Time|Round+Treatment+Bio_Rep, IRLcul_mumax, h = 7, quota = 1)
par(mfrow = c(3, 3))
par(mar = c(2, 2, 2, 1))
plot(IRLcul1 , log = "y")
summary(IRLcul1 )## $`1:Replete:1`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## -0.03861 0.16034 -0.15112 0.16012 -0.18606 -0.07377 0.12911
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.7592193 0.3438351 5.116 0.00372 **
## x 0.0075785 0.0006407 11.828 7.6e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1627 on 5 degrees of freedom
## Multiple R-squared: 0.9655, Adjusted R-squared: 0.9586
## F-statistic: 139.9 on 1 and 5 DF, p-value: 7.603e-05
##
##
## $`1:Replete:2`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## -0.05977 0.09664 0.02457 0.12515 -0.33697 0.07566 0.07473
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.793494 0.469100 1.692 0.151524
## x 0.007143 0.000691 10.336 0.000146 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1755 on 5 degrees of freedom
## Multiple R-squared: 0.9553, Adjusted R-squared: 0.9463
## F-statistic: 106.8 on 1 and 5 DF, p-value: 0.0001459
##
##
## $`1:Replete:3`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## -0.03760 -0.17761 0.53725 -0.39227 0.01679 0.10783 -0.05440
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.488243 0.895579 1.662 0.15745
## x 0.005008 0.001233 4.062 0.00971 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3132 on 5 degrees of freedom
## Multiple R-squared: 0.7674, Adjusted R-squared: 0.7209
## F-statistic: 16.5 on 1 and 5 DF, p-value: 0.009713
## y0 y0_lm mumax lag
## 1:Replete:1 74.22 5.807901 0.007578469 336.1912
## 1:Replete:2 88.61 2.211109 0.007142610 516.7230
## 1:Replete:3 71.52 4.429309 0.005007900 555.4690
## 1:Replete:1.r2 1:Replete:2.r2 1:Replete:3.r2
## 0.9654963 0.9552895 0.7674161
IRL culture trial 2 Replete:
IRLcul_mumax2 <- IRLcul_clean %>%
filter(Treatment %in% c("Replete"),
Round %in% c(2),
Time >= 0) %>%
rename(Chl = `Chlorophyll-A (RFU)`) %>%
arrange(Time, Bio_Rep)
IRLcul2 <- all_easylinear(Chl ~ Time|Round+Treatment+Bio_Rep, IRLcul_mumax2, h =7, quota = 1)
par(mfrow = c(3, 3))
par(mar = c(2, 2, 2, 1))
plot(IRLcul2 , log = "y")
summary(IRLcul2 )## $`2:Replete:1`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## -0.23499 0.44031 -0.16656 -0.01381 -0.09558 0.10722 -0.03659
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.7712453 0.4712913 5.880 0.002020 **
## x 0.0069299 0.0009628 7.198 0.000806 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2445 on 5 degrees of freedom
## Multiple R-squared: 0.912, Adjusted R-squared: 0.8944
## F-statistic: 51.81 on 1 and 5 DF, p-value: 0.0008063
##
##
## $`2:Replete:2`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## 0.15234 -0.55703 0.04284 0.55888 0.12964 -0.23617 -0.09050
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.626180 0.806476 2.016 0.09983 .
## x 0.008378 0.001503 5.575 0.00256 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3817 on 5 degrees of freedom
## Multiple R-squared: 0.8614, Adjusted R-squared: 0.8337
## F-statistic: 31.08 on 1 and 5 DF, p-value: 0.002558
##
##
## $`2:Replete:3`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## 0.31177 -0.51835 -0.13281 0.31474 0.18935 -0.07054 -0.09416
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.906698 0.628814 4.623 0.00572 **
## x 0.006738 0.001285 5.245 0.00334 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3263 on 5 degrees of freedom
## Multiple R-squared: 0.8462, Adjusted R-squared: 0.8155
## F-statistic: 27.51 on 1 and 5 DF, p-value: 0.003339
## y0 y0_lm mumax lag
## 2:Replete:1 73.64 15.978520 0.006929859 220.4869
## 2:Replete:2 81.69 5.084417 0.008378133 331.4284
## 2:Replete:3 71.42 18.296276 0.006738242 202.1121
## 2:Replete:1.r2 2:Replete:2.r2 2:Replete:3.r2
## 0.9119822 0.8614255 0.8462232
Remove outliers from OTB water data:
# Compute T0 baseline per replicate, treatment, and round
baseline <- aggregate(
`Chlorophyll-A (RFU)` ~ Bio_Rep + Treatment + Round,
data = subset(OTB, Time == 0),
FUN = mean,
na.rm = TRUE
)
names(baseline)[names(baseline) == "Chlorophyll-A (RFU)"] <- "Chl_T0"
# Merge baseline back into the full dataset
OTBwat_merged <- merge(OTB, baseline, by = c("Bio_Rep", "Treatment", "Round"), all.x = TRUE)
# Keep all T0 points, and remove any post-T0 points where replicate < its T0
OTBwat_clean <- subset(
OTBwat_merged,
Time == 0 | `Chlorophyll-A (RFU)` >= Chl_T0
)OTB water trial 1 Replete:
OTB_mumax1 <- OTBwat_clean %>%
filter(Treatment %in% c("Replete"),
Round %in% c(1),
Time >= 0) %>%
rename(Chl = `Chlorophyll-A (RFU)`) %>%
arrange(Time, Bio_Rep)
OTBwat1 <- all_easylinear(Chl ~ Time|Round+Treatment+Bio_Rep, OTB_mumax1, h = 7, quota = 1)
par(mfrow = c(3, 3))
par(mar = c(2, 2, 2, 1))
plot(OTBwat1 , log = "y")
summary(OTBwat1)## $`1:Replete:1`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## -0.15859 0.22093 -0.08966 -0.03358 0.16311 -0.01998 -0.08224
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.332309 0.210366 20.59 5.00e-06 ***
## x 0.007191 0.000602 11.95 7.25e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1529 on 5 degrees of freedom
## Multiple R-squared: 0.9661, Adjusted R-squared: 0.9594
## F-statistic: 142.7 on 1 and 5 DF, p-value: 7.251e-05
##
##
## $`1:Replete:2`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## 0.33247 -0.31312 -0.07991 -0.10915 0.06674 0.08438 0.01858
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.6685656 0.3405662 10.772 0.000120 ***
## x 0.0084842 0.0008604 9.861 0.000183 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2185 on 5 degrees of freedom
## Multiple R-squared: 0.9511, Adjusted R-squared: 0.9413
## F-statistic: 97.23 on 1 and 5 DF, p-value: 0.0001828
##
##
## $`1:Replete:3`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6 7
## 0.05109 0.05233 -0.25963 0.01366 0.20820 0.01297 -0.07862
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.3310777 0.2440626 17.75 1.04e-05 ***
## x 0.0068159 0.0006166 11.05 0.000106 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1566 on 5 degrees of freedom
## Multiple R-squared: 0.9607, Adjusted R-squared: 0.9528
## F-statistic: 122.2 on 1 and 5 DF, p-value: 0.0001055
## y0 y0_lm mumax lag
## 1:Replete:1 211.13 76.11986 0.007190613 141.8745
## 1:Replete:2 101.61 39.19564 0.008484247 112.2759
## 1:Replete:3 328.76 76.02618 0.006815913 214.8282
## 1:Replete:1.r2 1:Replete:2.r2 1:Replete:3.r2
## 0.9661413 0.9510923 0.9606888
OTB water trial 2 Replete:
OTB_mumax2 <- OTBwat_clean %>%
filter(Treatment %in% c("Replete"),
Round %in% c(2),
Time >= 0) %>%
rename(Chl = `Chlorophyll-A (RFU)`) %>%
arrange(Time, Bio_Rep)
OTBwat2 <- all_easylinear(Chl ~ Time|Round+Treatment+Bio_Rep, OTB_mumax2, h = 5, quota = 1)
par(mfrow = c(3, 3))
par(mar = c(2, 2, 2, 1))
plot(OTBwat2 , log = "y")
summary(OTBwat2)## $`2:Replete:1`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5
## 0.04205 -0.14349 0.12612 0.01002 -0.03470
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 5.0096618 0.0889572 56.315 1.23e-05 ***
## x 0.0041960 0.0007566 5.546 0.0116 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1148 on 3 degrees of freedom
## Multiple R-squared: 0.9111, Adjusted R-squared: 0.8815
## F-statistic: 30.76 on 1 and 3 DF, p-value: 0.01156
##
##
## $`2:Replete:2`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5
## 0.12058 0.06702 -0.35136 0.01932 0.14443
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.724234 0.180952 26.108 0.000123 ***
## x 0.007660 0.001539 4.977 0.015586 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2336 on 3 degrees of freedom
## Multiple R-squared: 0.892, Adjusted R-squared: 0.856
## F-statistic: 24.77 on 1 and 3 DF, p-value: 0.01559
##
##
## $`2:Replete:3`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5
## 0.1065 0.1344 -0.4463 0.0634 0.1420
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.409273 0.224859 19.609 0.00029 ***
## x 0.007442 0.001912 3.891 0.03009 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2903 on 3 degrees of freedom
## Multiple R-squared: 0.8347, Adjusted R-squared: 0.7795
## F-statistic: 15.14 on 1 and 3 DF, p-value: 0.03009
## y0 y0_lm mumax lag
## 2:Replete:1 156.29 149.85405 0.004195997 10.02180
## 2:Replete:2 127.08 112.64419 0.007660107 15.74165
## 2:Replete:3 91.45 82.20965 0.007442250 14.31284
## 2:Replete:1.r2 2:Replete:2.r2 2:Replete:3.r2
## 0.9111292 0.8919806 0.8346501
Remove outliers from IRL water data:
#Compute T0 baseline per replicate, treatment, and round
baseline <- aggregate(
`Chlorophyll-A (RFU)` ~ Bio_Rep + Treatment + Round,
data = subset(IRLrec, Time == 0),
FUN = mean,
na.rm = TRUE
)
names(baseline)[names(baseline) == "Chlorophyll-A (RFU)"] <- "Chl_T0"
# Merge baseline back into the full dataset
IRLwat_merged <- merge(IRLrec, baseline, by = c("Bio_Rep", "Treatment", "Round"), all.x = TRUE)
# Keep all T0 points, and remove any post-T0 points where replicate < its T0
IRLwat_clean <- subset(
IRLwat_merged,
Time == 0 | `Chlorophyll-A (RFU)` >= Chl_T0
)IRL water trial 1 Replete:
IRL_mumax1 <- IRLwat_clean %>%
filter(Treatment %in% c("Replete"), Round %in% c(1)) %>%
rename(Chl = `Chlorophyll-A (RFU)`) %>%
droplevels()
IRLwat1 <- all_easylinear(Chl ~ Time|Round+Treatment+Bio_Rep, IRL_mumax1, h = 5, quota = 1)
par(mfrow = c(3, 3))
par(mar = c(2, 2, 2, 1))
plot(IRLwat1 , log = "y")
summary(IRLwat1)## $`1:Replete:1`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5
## -0.07587 0.09142 -0.13804 -0.22946 0.35195
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 5.005693 0.268748 18.626 0.000338 ***
## x 0.008304 0.001006 8.258 0.003718 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2644 on 3 degrees of freedom
## Multiple R-squared: 0.9579, Adjusted R-squared: 0.9438
## F-statistic: 68.2 on 1 and 3 DF, p-value: 0.003718
##
##
## $`1:Replete:2`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5
## -0.01569 -0.18763 0.22688 -0.29711 0.27355
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 5.6061505 0.2231989 25.12 0.000138 ***
## x 0.0064232 0.0008246 7.79 0.004403 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2887 on 3 degrees of freedom
## Multiple R-squared: 0.9529, Adjusted R-squared: 0.9372
## F-statistic: 60.68 on 1 and 3 DF, p-value: 0.004403
##
##
## $`1:Replete:3`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5
## -0.32328 0.13706 0.35205 -0.08458 -0.08125
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.230102 0.228469 18.515 0.000344 ***
## x 0.013307 0.001943 6.848 0.006374 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.295 on 3 degrees of freedom
## Multiple R-squared: 0.9399, Adjusted R-squared: 0.9198
## F-statistic: 46.89 on 1 and 3 DF, p-value: 0.006374
## y0 y0_lm mumax lag
## 1:Replete:1 108.14 149.26049 0.008304417 -38.806618
## 1:Replete:2 267.86 272.09478 0.006423229 -2.442075
## 1:Replete:3 178.44 68.72424 0.013306580 71.705157
## 1:Replete:1.r2 1:Replete:2.r2 1:Replete:3.r2
## 0.9578639 0.9528881 0.9398718
IRL water trial 2 Replete:
IRL_mumax2 <- IRLwat_clean %>%
filter(Treatment %in% c("Replete"),
Round %in% c(2),
Time >= 0) %>% # <-- exclude T0
rename(Chl = `Chlorophyll-A (RFU)`) %>%
arrange(Time, Bio_Rep)
IRLwat2 <- all_easylinear(Chl ~ Time|Round+Treatment+Bio_Rep, IRL_mumax2, h =6, quota = 1)
par(mfrow = c(3, 3))
par(mar = c(2, 2, 2, 1))
plot(IRLwat2 , log = "y")
summary(IRLwat2)## $`2:Replete:1`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6
## 0.27656 -0.39465 -0.12399 0.08819 0.39137 -0.23748
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.126211 0.317395 13.000 0.000202 ***
## x 0.009866 0.001698 5.811 0.004365 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3409 on 4 degrees of freedom
## Multiple R-squared: 0.8941, Adjusted R-squared: 0.8676
## F-statistic: 33.76 on 1 and 4 DF, p-value: 0.004365
##
##
## $`2:Replete:2`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6
## 0.03290 0.04031 -0.06486 0.15651 -0.44416 0.27930
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.883134 0.200405 24.366 1.68e-05 ***
## x 0.006915 0.001379 5.014 0.00741 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2769 on 4 degrees of freedom
## Multiple R-squared: 0.8628, Adjusted R-squared: 0.8284
## F-statistic: 25.14 on 1 and 4 DF, p-value: 0.007414
##
##
## $`2:Replete:3`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5 6
## -0.10408 0.24742 0.10827 -0.59427 0.39445 -0.05179
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.516376 0.359111 12.577 0.00023 ***
## x 0.008128 0.001921 4.231 0.01336 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3857 on 4 degrees of freedom
## Multiple R-squared: 0.8174, Adjusted R-squared: 0.7717
## F-statistic: 17.9 on 1 and 4 DF, p-value: 0.01336
## y0 y0_lm mumax lag
## 2:Replete:1 102.18 61.94276 0.009865993 50.732380
## 2:Replete:2 136.46 132.04390 0.006914858 4.757449
## 2:Replete:3 105.76 91.50335 0.008128438 17.813607
## 2:Replete:1.r2 2:Replete:2.r2 2:Replete:3.r2
## 0.8940791 0.8627528 0.8173779
IRL water trial 1 -B12:
IRL_mumax1_b <- IRLwat_clean %>%
filter(Treatment %in% c("-Cobalamin"),
Round %in% c(1),
Time >= 0) %>% # <-- exclude T0
rename(Chl = `Chlorophyll-A (RFU)`) %>%
arrange(Time, Bio_Rep)
IRLwat1_b <- all_easylinear(Chl ~ Time|Round+Treatment+Bio_Rep, IRL_mumax1_b, h = 4, quota =1)
par(mfrow = c(3, 3))
par(mar = c(2, 2, 2, 1))
plot(IRLwat1_b , log = "y")
summary(IRLwat1_b)## $`1:-Cobalamin:1`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4
## -0.13839 0.02356 0.36806 -0.25322
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.629046 0.277043 16.709 0.00356 **
## x 0.011942 0.003085 3.871 0.06073 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3311 on 2 degrees of freedom
## Multiple R-squared: 0.8822, Adjusted R-squared: 0.8233
## F-statistic: 14.98 on 1 and 2 DF, p-value: 0.06073
##
##
## $`1:-Cobalamin:2`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4
## 0.09071 -0.12902 -0.06650 0.10480
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 5.1181197 0.1239394 41.295 0.000586 ***
## x 0.0053753 0.0007698 6.982 0.019900 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1419 on 2 degrees of freedom
## Multiple R-squared: 0.9606, Adjusted R-squared: 0.9409
## F-statistic: 48.76 on 1 and 2 DF, p-value: 0.0199
##
##
## $`1:-Cobalamin:3`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4
## 0.06389 -0.33062 0.46957 -0.20284
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.332327 0.711542 4.683 0.0427 *
## x 0.015666 0.004035 3.883 0.0604 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.433 on 2 degrees of freedom
## Multiple R-squared: 0.8829, Adjusted R-squared: 0.8243
## F-statistic: 15.08 on 1 and 2 DF, p-value: 0.06038
## y0 y0_lm mumax lag
## 1:-Cobalamin:1 89.18 102.41631 0.011941636 -11.58880
## 1:-Cobalamin:2 182.88 167.02102 0.005375287 16.87548
## 1:-Cobalamin:3 112.12 28.00342 0.015665946 88.55150
## 1:-Cobalamin:1.r2 1:-Cobalamin:2.r2 1:-Cobalamin:3.r2
## 0.8822319 0.9605951 0.8828865
IRL water trial 2 -B12:
IRL_mumax2_b <- IRLwat_clean %>%
filter(Treatment %in% c("-Cobalamin"),
Round %in% c(2),
Time >= 0) %>%
rename(Chl = `Chlorophyll-A (RFU)`) %>%
arrange(Time, Bio_Rep)
IRLwat2_b <- all_easylinear(Chl ~ Time|Round+Treatment+Bio_Rep, IRL_mumax2_b, h =5, quota =1)
par(mfrow = c(3, 3))
par(mar = c(2, 2, 2, 1))
plot(IRLwat2_b , log = "y")
summary(IRLwat2_b)## $`2:-Cobalamin:1`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5
## 0.08915 -0.35337 0.17809 0.34733 -0.26120
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.371518 0.460476 9.493 0.00248 **
## x 0.008698 0.002261 3.847 0.03101 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3432 on 3 degrees of freedom
## Multiple R-squared: 0.8314, Adjusted R-squared: 0.7752
## F-statistic: 14.8 on 1 and 3 DF, p-value: 0.03101
##
##
## $`2:-Cobalamin:2`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5
## -0.04172 0.31220 -0.17814 -0.37176 0.27941
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 5.022854 0.290615 17.284 0.000422 ***
## x 0.007184 0.001842 3.899 0.029934 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3402 on 3 degrees of freedom
## Multiple R-squared: 0.8352, Adjusted R-squared: 0.7803
## F-statistic: 15.2 on 1 and 3 DF, p-value: 0.02993
##
##
## $`2:-Cobalamin:3`
##
## Call:
## lm(formula = y ~ x)
##
## Residuals:
## 1 2 3 4 5
## 0.23747 0.05598 -0.51117 -0.09547 0.31319
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.206086 0.396161 10.617 0.00179 **
## x 0.010498 0.002488 4.219 0.02434 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3777 on 3 degrees of freedom
## Multiple R-squared: 0.8558, Adjusted R-squared: 0.8077
## F-statistic: 17.8 on 1 and 3 DF, p-value: 0.02434
## y0 y0_lm mumax lag
## 2:-Cobalamin:1 107.43 79.16373 0.008697853 35.103053
## 2:-Cobalamin:2 145.64 151.84401 0.007183847 -5.806909
## 2:-Cobalamin:3 130.09 67.09340 0.010498480 63.070164
## 2:-Cobalamin:1.r2 2:-Cobalamin:2.r2 2:-Cobalamin:3.r2
## 0.8314296 0.8352092 0.8557594
Extract and combine all specific growth rates:
# Helper function to handle matrix or vector coefficients
make_df <- function(model, round_num, experiment_name) {
coefs <- coef(model)
if (is.matrix(coefs)) {
df <- as.data.frame(coefs) %>%
tibble::rownames_to_column("Group") %>%
mutate(
Round = round_num,
Experiment = experiment_name,
r2 = rsquared(model)
)
} else {
df <- data.frame(
Group = names(coefs),
Estimate = as.numeric(coefs),
Round = round_num,
Experiment = experiment_name,
r2 = rep(rsquared(model), length(coefs))
)
}
return(df)
}
# Create data frames for all models with experiment names
df1 <- make_df(IRLcul1, 1, "IRL Culture")
df2 <- make_df(IRLcul2, 2, "IRL Culture")
df3 <- make_df(OTBwat1, 1, "OTB Water")
df4 <- make_df(OTBwat2, 2, "OTB Water")
df5 <- make_df(IRLwat1, 1, "IRL Water")
df6 <- make_df(IRLwat2, 2, "IRL Water")
df7 <- make_df(IRLwat1_b, 1, "IRL Water B")
df8 <- make_df(IRLwat2_b, 2, "IRL Water B")
# Combine all results
all_results <- bind_rows(df1, df2, df3, df4, df5, df6, df7, df8)
# Optional: reorder columns
all_results <- all_results %>%
select(Experiment, Round, Group, y0, mumax, r2)Prepare IRL seawater data for T-test for both trials in -B12 media:
df5 <- as.data.frame(coef(IRLwat1))
df5$Round <- 1
df5$Experiment <- "IRL Water"
df6 <- as.data.frame(coef(IRLwat2))
df6$Round <- 2
df6$Experiment <- "IRL Water"
df7 <- as.data.frame(coef(IRLwat1_b))
df7$Round <- 1
df7$Experiment <- "IRL Water -B12"
df8 <- as.data.frame(coef(IRLwat2_b))
df8$Round <- 2
df8$Experiment <- "IRL Water -B12"
# Combine
IRLwater <- bind_rows(df5, df6, df7, df8)T-test between specific growthrates both trials for IRL water in -B12:
##
## Two Sample t-test
##
## data: mumax by Experiment
## t = -0.59074, df = 10, p-value = 0.5678
## alternative hypothesis: true difference in means between group IRL Water and group IRL Water -B12 is not equal to 0
## 95 percent confidence interval:
## -0.005105451 0.002965605
## sample estimates:
## mean in group IRL Water mean in group IRL Water -B12
## 0.008823919 0.009893842
Plot IRL seawater specific growth rates:
# Extract the columns needed
mumax_results <- all_results %>%
select(Experiment, Round, Group, mumax) %>%
mutate(ExperimentRound = paste(Experiment, Round, sep = "_"))
mumax_results <- mumax_results %>%
mutate(
Treatment = factor(
Experiment,
levels = c("IRL Water B", "IRL Water") # order: B12 first, Replete second
)
)
# Filter IRL Water only
mumax_IRLwater <- mumax_results %>%
filter(grepl("IRL Water", Experiment))
IRL_mu <- ggplot(mumax_IRLwater, aes(x = Treatment, y = mumax, color = Treatment)) +
# Raw replicate points
geom_point(size = 5, shape = 16, position = position_jitter(width = 0.1)) +
# Mean per Experiment+Round (open circle)
stat_summary(fun = mean, geom = "point", shape = 1, size = 5, stroke = 1.2) +
stat_summary(fun.data = mean_se, geom = "errorbar", width = 0.2) +
# facet_grid(~Round) +
scale_color_manual(values = c(
"IRL Water B" = "#102E50",
"IRL Water" = "#F5C45E",
"-Thiamine" = "deeppink4",
"-Biotin" = "#667028"
)) +
scale_x_discrete(labels = c(
"IRL Water B" = expression("\u2013B"[12]),
"IRL Water" = "Replete",
"-Thiamine" = expression("\u2013Thiamine"),
"-Biotin" = expression("\u2013Biotin")
))+
theme_minimal() +
theme(
legend.position = "none",
panel.grid = element_blank(),
panel.border = element_rect(fill = NA, color = "black"),
axis.title = element_text(size = 20),
axis.text = element_text(size = 20, color = "black"),
axis.text.x = element_text( size = 18, color = "black", angle = 45, hjust = 1),
strip.text = element_text( size = 20, color = "black"),
legend.title = element_blank(),
legend.text = element_blank(),
legend.key = element_blank(),
axis.ticks = element_line(color = "black", linewidth = 0.5),
axis.ticks.length = unit(4, "pt")
) +
labs(
x = "Media Treatment",
y = expression(paste("(", mu, ") RFU" * "\u00B7" * "day"^-1, " "))
) +
guides(color = "none", fill = "none", shape = "none") +
# Set visible range safely
coord_cartesian(ylim = c(0, 0.5))
IRL_mu# Ensure BOTH subplots have no legend inside patchwork
max_bio__2 <- max_bio_plot + theme(legend.position = "none")
mu_2 <- IRL_mu + theme(legend.position = "none")
bo <-
growth /
(max_bio_plot +
IRL_mu ) +
plot_layout(
guides = "collect",
widths = c(4, 2)
) &
theme(
legend.position = "bottom",
panel.spacing = unit(.1, "lines"),
plot.margin = margin(15, 15, 15, 15),
axis.title.y = element_text(margin = margin(r = 12)),
axis.title.x = element_text(margin = margin(t = 8)),
legend.margin = margin(t = 4, b = 4),
legend.box.margin = margin(5, 5, 5, 5),
axis.text.x = element_text(size = 16)
)
boImport data:
# Identify sample columns
sample_cols <- names(Recruit_16s)[2:(which(names(Recruit_16s) == "total") - 1)]
# Pivot to long data
Recruit_16s_long <- Recruit_16s %>%
pivot_longer(
cols = all_of(sample_cols),
names_to = "sample",
values_to = "count"
)
# Separate samples into new columns
Recruit_16s_long <- Recruit_16s_long %>%
separate(sample, into = c("experiment","treatment","filter_size","bio_rep"),
sep = "_", remove = FALSE)
# One genus in one sample with the total count of that genus in that sample
grouped_data1 <- Recruit_16s_long %>%
mutate(count = as.numeric(count)) %>%
group_by(treatment, experiment, filter_size, bio_rep, genus) %>%
summarise(total_count = sum(count), .groups = "drop"
)Prepare data as wide matrix
# Pivot wider
wide_data1 <- grouped_data1 %>%
pivot_wider(names_from = genus, values_from = total_count, values_fill = 0)
# Create rownames and remove metadata columns for distance matrix
rownames(wide_data1) <- paste(wide_data1$treatment, wide_data1$experiment, wide_data1$filter_size, wide_data1$bio_rep, sep = "_")
wide_matrix1 <- wide_data1[, !(names(wide_data1) %in% c("treatment", "experiment", "filter_size", "bio_rep"))]
wide_matrix1 <- as.data.frame(lapply(wide_matrix1, as.numeric))
wide_matrix1[is.na(wide_matrix1)] <- 0
# Bray-Curtis + PCoA
bray_dist1 <- vegdist(wide_matrix1, method = "bray")
pcoa_result1 <- pcoa(bray_dist1)
# Build plot data
pcoa_data1 <- as.data.frame(pcoa_result1$vectors)
pcoa_data1$treatment <- wide_data1$treatment
pcoa_data1$experiment <- wide_data1$experiment
pcoa_data1$filter_size <- as.factor(wide_data1$filter_size)Permanova:
#ensure wide matrix is numeric
wide_matrix1 <- as.data.frame(lapply(wide_matrix1, as.numeric))
wide_matrix1[is.na(wide_matrix1)] <- 0
#Build the Bray–Curtis distance matrix
bray_dist <- vegdist(wide_matrix1, method = "bray")
#Bring metadata back in (from wide_data1)
metadata <- wide_data1[, c("treatment", "experiment", "filter_size", "bio_rep")]
permanova_result <- adonis2(
bray_dist ~ experiment,
data = metadata,
permutations = 999
)
print(permanova_result)## Permutation test for adonis under reduced model
## Permutation: free
## Number of permutations: 999
##
## adonis2(formula = bray_dist ~ experiment, data = metadata, permutations = 999)
## Df SumOfSqs R2 F Pr(>F)
## Model 4 6.2820 0.67542 16.127 0.001 ***
## Residual 31 3.0189 0.32458
## Total 35 9.3009 1.00000
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
sw_pcoa <- ggplot(
pcoa_data1,
aes(
x = Axis.1,
y = Axis.2,
color = experiment,
shape = treatment,
alpha = filter_size # keep as factor
)
) +
geom_point(size = 5) +
labs(
x = paste0("PCoA1 (", round(pcoa_result1$values$Relative_eig[1] * 100, 2), "%)"),
y = paste0("PCoA2 (", round(pcoa_result1$values$Relative_eig[2] * 100, 2), "%)"),
color = "Experiment",
shape = "Treatment",
alpha = "Filter size"
) +
stat_ellipse(aes(group = experiment)) +
scale_color_hue() +
# Map discrete alpha manually
scale_alpha_manual(values = c("02" = 0.4, "10" = 1)) +
theme_minimal() +
theme(
legend.key.size = unit(1.5, "cm"),
legend.title = element_text(size = 20),
legend.text = element_text(size = 20),
panel.border = element_rect(color = 'black', fill = NA),
axis.title.x = element_text(size = 20, color = "black"),
axis.title.y = element_text(size = 20, color = "black"),
strip.text = element_text(size = 20, color = "black"),
axis.text.y = element_text(size = 20, color = "black"),
axis.text.x = element_text(size = 20, color = "black")
) +
geom_hline(yintercept = 0, linetype = "dotted", linewidth = 0.5) +
geom_vline(xintercept = 0, linetype = "dotted", linewidth = 0.5) +
theme(panel.grid = element_blank())
plot(sw_pcoa)Prepare data for relative abundance plot:
#Standardize experiment names
grouped_data <- grouped_data1 %>%
mutate(
experiment = case_when(
grepl("1003.1006", experiment, ignore.case = TRUE) ~ "1003.1006",
grepl("OTB", experiment, ignore.case = TRUE) ~ "OTB.1006",
grepl("IRL", experiment, ignore.case = TRUE) ~ "IRL.1006",
TRUE ~ experiment
),
experiment = factor(experiment,
levels = c("1003.1006", "OTB_seawater", "OTB.1006", "IRL_seawater", "IRL.1006"),
labels = c("Microbiome from IRL culture",
"OTB Seawater",
"Microbiome from OTB seawater",
"IRL seawater",
"Microbiome from IRL seawater"))
)
#Pivot to long format
sample_cols <- names(Recruit_16s)[2:(which(names(Recruit_16s) == "total") - 1)]
Recruit_16s_long <- Recruit_16s %>%
pivot_longer(
cols = all_of(sample_cols),
names_to = "sample",
values_to = "count"
) %>%
separate(sample, into = c("experiment","treatment","filter_size","bio_rep"),
sep = "_", remove = FALSE) %>%
mutate(count = as.numeric(count)) # ensure numeric
# ️ Summarize counts per sample & genus
grouped_data <- Recruit_16s_long %>%
group_by(experiment, treatment, filter_size, bio_rep, order) %>%
summarise(total_count = sum(count), .groups = "drop") %>%
group_by(experiment, treatment, filter_size, bio_rep) %>%
mutate(RA = (total_count / sum(total_count)) * 100) %>%
ungroup() %>%
filter(treatment %in% c("noB12", "replete", "seawater"))
# ⃣Collapse low abundance genera (<1%) into "z<1% abund."
grouped_data$order <- as.character(grouped_data$order)
grouped_data$order [grouped_data$RA < 1] <- "z<1% abund."
grouped_data$order <- factor(grouped_data$order ) # back to factor
# Optional: fix experiment & filter labels for plotting
grouped_data <- grouped_data %>%
mutate(
experiment = factor(experiment,
levels = c("X1003.1006","OTB","OTB.1006","IRL","IRL.1006"),
labels = c("Recruited IRL culture",
"OTB Seawater",
"Recruited OTB seawater",
"IRL Seawater",
"Recruited IRL seawater")),
Treatment = recode(treatment, "-B12" = "\u2014B12"),
filter_size = recode(filter_size, "10" = "10.0um", "02" = "0.2um"),
filter_size = factor(filter_size, levels = c("10.0um", "0.2um")), # <-- top first
x_label = paste(Treatment, "Rep", bio_rep, sep = " ")
)
# Define color palette (ensure enough colors for all genera)
unique_genera <- levels(grouped_data$order )
# Use your custom palette if available, otherwise RColorBrewer
library(RColorBrewer)
n_colors <- length(unique_genera)
palette_colors <- brewer.pal(min(12, n_colors), "Set3") # adjust if >12 genera
if(n_colors > 12) {
palette_colors <- colorRampPalette(palette_colors)(n_colors)
}
names(palette_colors) <- unique_generataxa_bar_plot <- ggplot(grouped_data,
aes(x = x_label, y = RA, fill = order)) +
geom_bar(stat = "identity", position = "stack") +
facet_grid(filter_size ~ experiment, scales = "free_x") +
theme_classic() +
theme(
legend.position = "right",
legend.title = element_blank(),
axis.text.x = element_text(size = 18, angle = 45, hjust = 1, color = "black"),
axis.text.y = element_text(size = 18, color = "black"),
axis.title.y = element_text(size = 18),
strip.text = element_text(size = 12),
legend.text = element_text(size = 18),
strip.position = "top"
) +
scale_y_continuous() + # Corrected
scale_fill_viridis_d(option = "turbo") + # Deep distinct colors
guides(fill = guide_legend(ncol = 2)) +
ylab("Relative Abundance (%)") +
labs(x = NULL)
taxa_bar_plot